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1.  VLP  RADIO  HATE  PROPAGATION  TO  GREAT  DISTMKXS 

SCOPE  OF  THIS  DOCUMENT 

The  Omega  Navigation  System  Operations  Detail  (ONSOD)  of  the  US  Coast 
Guard  is  responsible  for  the  implementation  and  evaluation  of  the  Omega 
Navigation  System.  Omega  uses  very  low  frequency  (VLF)  radio  waves  with  range 
from  10.2  to  13.6  kHz  and  is  an  operational  global  system.  Thus,  a  complete 
evaluation  of  the  system  requires  an  understanding  of  VLF  radio  wave  propaga¬ 
tion  under  both  spatial  and  tenqporal  operating  conditions.  A  particular  set 
of  conditions  has  been  of  considerable  interest,  namely  transequatorial  propa¬ 
gation  at  night  or  during  day-night  transitions.  Theoretical  evaluation  of 
such  propagation  paths  requires  highly  detailed  tracking  of  the  propagating 
modes  through  the  region  of  the  geomagnetic  equator  because  the  modes  can  be 
very  sensitive  to  the  magnetic  field.  An  even  more  difficult  problem  is  the 
analysis  of  the  effect  of  a  moving  day-night  terminator  on  such  paths. 

A  general  background  review  of  VLF  radio  wave  propagation  from  a  theo¬ 
retical  viewpoint  is  presented,  and  this  is  followed  by  a  discussion  of  trans¬ 
equatorial  VLF  propagation  at  night  and  during  day-night  transitions.  A  VLF 
propagation  model  which  includes  consideration  of  inhomogeneities  both  in  the 
vertical  direction  and  along  the  propagation  path  is  then  presented.  The 
model  includes  consideration  of  the  anisotropic  effects  (caused  by  the  geomag¬ 
netic  field) .  It  has  been  implemented  into  a  FORTRAN  coiqputer  program  called 
FASTMC.  This  program  allows  for  an  arbitrarily  oriented  and  elevated  trans¬ 
mitting  dipole  antenna  from  which  the  vertical  or  horizontal  electric  field  at 
an  arbitrary  elevation  can  be  output.  An  application  of  the  FASTMC  capabili¬ 
ties  has  been  further  i^>lemented  into  a  program  called  BUMP,  which  allows  for 
relatively  easy  modeling  of  a  moving  perturbation  with  a  minimum  of  data  prep¬ 
aration.  The  FASTMC- BUMP  combination  is  designed  specifically  for  consider¬ 
ation  of  transequatorial  VLF  propagation.  A  discussion  of  both  FASTMC  and 
BUMP  is  included. 

THEORETICAL  DEVELOPMENTS 

VLF  propagation  to  great  distances  can  be  conveniently  represented  in 
terms  of  waveguide  node  propagation,  where  the  finitely  conducting  curved 
earth  and  the  anisotropic,  imperfectly  conducting  curved  ionosphere  with  dip- 


ping  Magnetic  field  form  the  boundaries  of  a  waveguide.  The  basic  idea  that 
the  earth  and  lower  ionosphere  form  a  waveguide  actually  goes  back  to  the  work 
of  Watson  (1918,  1919),  who  postulated  that  the  Kennelly-Heaviside  layer  could 
be  represented  as  a  conducting  shell  concentric  with  the  spherical  conducting 
earth.  Watson  developed  an  exact,  slowly  convergent,  harmonic  series-type 
solution  to  the  EM  propagation  problem  and  further  devised  a  transformation  of 
the  series  to  a  highly  convergent  form.  The  transformation,  which  involves 
distortion  of  the  contour  of  integration  in  the  complex  wave  number  plane,  is 
now  identified  with  Watson' s  name,  although  the  transformed  series  is  often 
"rediscovered"  by  subsequent  workers.  Although  "waveguides"  were  not  really 
"discovered"  until  the  1930s,  the  terms  in  Watson's  transformed  series  were 
simply  waveguide  modes. 

The  theory  of  VLF  waveguide  mode  propagation  has  been  presented  in  the 
literature  in  two  different  formulations.  One  formulation,  based  on  the 
pioneering  work  of  Watson  (1919)  and  more  recently  in  the  theory  given  by 
Schumann  (1952,  1954)  and  developed  extensively  by  Wait  (1962),  treats  the 
problem  in  spherical  coordinates  wherein  the  fields  are  expanded  in  terms  of 
azimuthal,  longitudinal,  and  radial  functions.  In  a  sequence  of  simplifying 
approximations,  including  the  assunption  of  azimuthal  and  longitudinal  inde¬ 
pendence  of  both  ground  and  ionospheric  properties.  Wait  (1962)  develops  a 
modal  (eigenvalue)  equation  in  terms  of  Airy  functions.  The  Airy  functions 
come  about  as  a  result  of  a  third-order  approximation  to  spherical  wave 
functions  (Watson,  1944) .  The  solution  of  the  eigenvalue  equation  yields  the 
propagation  characteristics  of  the  VLF  modes  in  the  earth-ionosphere  wave¬ 
guide,  and  consideration  of  orthogonality  properties  of  the  radial  (Airy) 
functions  yields  excitation  factors  of  the  modes.  The  excitation  factor  is  a 
quantity  that  gives  the  aaplitude  of  the  wave  excited  in  a  given  mode  by  a 
given  source.  Another  sinplifying  approximation  often  made  involves  the 
replacement  of  the  inhomogeneous  ionosphere  with  an  impedance  boundary  placed 
at  the  height  where  the  bulk  of  the  VLF  energy  is  assumed  to  be  reflected.  A 
similar  iapedance  boundary  replacement  is  also  made  for  the  earth.  The  impe¬ 
dance  boundary  formulation  permits  specification  of  the  ratio  of  an  electric 
field  cooponent  to  a  magnetic  field  conponent  at  the  boundary  without  further 
consideration  of  the  region  beyond  the  boundary.  Analytically,  formulation  in 
terms  of  Impedance  boundaries  is  equivalent  to  assumption  of  "homogeneous 
boundary  conditions"  (Friedman,  1959;  Morse  and  Feshbach,  1953).  The  replace- 
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ment  of  an  Inhomogeneous*  isotropic  ionosphere  by  an  "equivalent"  impedance 
boundary  is  a  straightforward  procedure.  Such  a  replacement  for  an  anisotrop¬ 
ic  ionosphere  is  greatly  conf>licated  by  the  coupling  of  transverse  magnetic 
and  transverse  electric  polarizations.  Boundary  conditions  oust  then  be 
formulated  in  terms  of  an  inpedance  matrix.  unfortunately,  the  conplication 
is  sufficient  to  render  many  of  Wait's  results  inapplicable  to  the  highly 
anisotropic  nighttime  ionosphere  and  to  limit  their  application  to  the  daytime 
ionosphere,  for  which  the  effect  of  anisotropy  is  slight.  The  theory  has  been 
significantly  modified  and  extended  to  anisotropic  ionospheres  by  Galejs 
(1972). 

Another  formulation  of  the  earth-ionosphere  waveguide  problem  has  been 
presented  by  Budden  (1961a).  In  this  formulation,  a  modal  equation  is  devel¬ 
oped  in  terms  of  reflection  coefficients  of  the  ionosphere  and  the  earth. 
Formulation  in  terms .  of  reflection  coefficients  permits  determination  of 
propagation  characteristics  of  waveguide  modes  for  essentially  any  ionospheric 
or  ground  conditions,  limited  only  by  the  ability  to  compute  reflection  coef¬ 
ficients  for  the  environment  considered.  Budden' s  formulation  is  essentially 
planar,  although  he  does  point  out  how  to  include  earth  curvature  in  the 
direction  of  propagation  in  an  approximate  way  by  modifying  the  refractive 
index  in  the  space  between  the  earth  and  the  ionosphere  (Budden,  1962).  A 
more  rigorous  technique  was  introduced  by  Richter  (1966)  and  has  been  employed 
in  VLF  propagation  studies  by  Pappert  (1968)  and  Abbas  et  al  (1971).  This 
technique  assumes  cylindrical  stratification  and  effects  a  conformal  transfor¬ 
mation  from  cylindrical  to  Cartesian  coordinates.  Modal  excitation  values  are 
determined  in  Budden' s  formulation  by  applying  of  refidue  theory  rather  than 
by  using  the  modal  orthogonality  properties  of  conventional  waveguide  theory 
employed  by  Wait. 

Whichever  formulation  is  used  at  the  outset,  the  use  of  waveguide  mode 
theory  leads  to  specification  of  the  VLF  field  as  a  summation  of  terms  called 
inodes.  Because  most  VLF  transmitters  in  common  use  radiate  a  vertically 
polarized  field,  only  the  radial  or  vertical  component  of  the  electric  field 
is  usually  considered.  For  a  time  harmonic  bource,  the  VLF  mode  Bum  for  the 
vertical  electric  field  may  be  written  as 
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where  K(P,f)  is  a  corn? lex  constant  dependent  on  transmitted  power  (P)  and 
frequency  (f ) ,  d  is  the  distance  from  the  transmitter  on  a  homogeneous,  smooth 
earth  of  radius  a,  AR  is  the  excitation  factor  for  mode  n,  normalized  to  unity 
for  flat  earth,  perfectly  conducting  boundaries,  kQ  is  the  free  space  wave 
number,  Sn  is  the  propagation  factor,  and  Gn'  represents  height-gain  func¬ 
tions  for  mode  n,  normalized  to  unity  at  the  ground.  One  height-gain  function 
is  needed  for  the  transmitter  (T)  and  one  for  the  receiver  (R).  Generally, 
both  AR  and  SR  are  complex.  The  real  part  of  Sn  determines  the  distance 
dependence  of  the  phase  for  a  mode  while  the  imaginary  part  of  SR  determines 
the  attenuation  rate.  Thus 
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where  is  the  modal  attenuation  rate,  vR  is  the  phase  velocity,  and  c  is  the 
speed  of  light  in  free  space. 

Numerous  assumptions  and  approximations,  some  of  which  will  be  removed  in 
following  sections,  are  implicit  in  equation  ( 1 ) .  It  is  assumed  that  the 
propagation  environment  is  homogeneous  in  all  but  the  vertical  direction. 
Thus,  the  propagation  factor,  SR,  is  independent  of  position.  The  influence 
of  the  curvature  of  the  earth  transverse  to  the  propagation  direction  is 
approximately  accounted  for  in  the  spreading  term,  sin(d/a).  The  energy  flux 
for  a  wave  traveling  outwards  with  cylindrical  symmetry  on  a  flat  earth  with¬ 
out  attenuation  is  proportional  to  1/r.  For  a  spherically  curved  earth  this 
energy  flux  would  be  proportional  to  1/sin(r/a).  Thus  the  usual  1//r  term  for 
a  cylindrical  wave  is  replaced  by  the  1//sin(r/a)  term.  A  careful  derivation 
of  the  mode  equation  in  spherical  coordinates,  such  as  that  presented  by  Wait 
(1962)  and  discussed  previously,  yields  the  sin(d/a)  term  as  a  consequence  of 
an  asymptotic  approximation  to  a  Legendre  function.  The  mode  sum  in  equation 
(1)  further  assumes  that  "round  the  world"  signals  are  absent,  an  approxima- 
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tion  not  valid  when  receiver  locations  are  near  the  antipode  of  the  transmit¬ 
ter  —  multipath  focusing  then  occurs. 

Extensive  results  of  numerical  calculations  of  VLF  waveguide  mode  con¬ 
stants  have  been  presented  in  the  literature.  Early  results  given  by  Wait  and 
Spies  (1964)  show  some  general  properties  of  VLF  modes  for  typical  daytime  and 
nighttime  ionospheres,  with  various  ground  conductivities  and  for  east-to-west 
or  west-to-east  propagation  directions  at  the  magnetic  equator  or  when  the 
earth's  magnetic  field  is  ignored.  In  general,  the  attenuation  rate  of  the 
lowest  order  mode  is  less  at  night  than  during  the  day,  and  is  less  for  west- 
to-east  propagation  than  for  east-to-west.  The  attenuation  rate  increases  as 
ground  conductivity  decreases  —  rapidly  at  low  conductivities  such  as  occur 
for  the  polar  icecaps.  Each  successive  mode  order  suffers  a  greater  attenua¬ 
tion  rate  than  the  previous  order  mode  for  both  day  and  night. 

Typically,  daytime  attenuation  rates  for  the  first  order  mode  range  from 
2  to  4  dB/1000  km  over  the  VLF  band  for  sea  water  ground  conditions  (essen¬ 
tially  equivalent  to  a  perfect  conductor  at  VLF).  The  nighttime  attenuation 
rates  range  from  about  0.5  to  2  dB/1000  km  over  the  VLF  band  for  the  same 
ground  conductivities. 

The  results  of  Wait  and  Spies  show  that  the  modal  phase  velocities  for 
the  curved  earth  can  be  less  than  the  free  space  speed  of  light  at  intermedi¬ 
ate  to  higher  VLF  frequencies  but  increase  markedly  at  lower  frequencies. 
Horeover,  phase  velocities  increase  at  low  ionosphere  heights  (daytime)  com¬ 
pared  to  hitler  heights  (night),  and  a  decrease  in  ground  conductivity 
produces  a  decrease  in  phase  velocity.  An  increase  in  phase  velocity  with 
mode  number  is  also  indicated. 

The  modal  excitation  factors  and  height  gains  given  by  Wait  and  Spies 
show  little  dependence  on  the  terrestrial  magnetic  field.  At  higher  frequen¬ 
cies  for  isotropic  conditions,  the  excitation  factor  for  the  first  mode  is 
significantly  reduced  in  comparison  with  that  for  a  flat  earth  case.  The 
second  mode  exhibits  a  ouch  smaller  reduction  at  higher  frequencies.  Gener¬ 
ally,  lower  ground  conductivities  produce  an  increase  in  excitation  factor. 
The  excitation  factor  for  the  second  order  mode  is  significantly  higher  than 
that  for  the  first  order  mode,  especially  at  higher  frequencies  and  at  night. 
At  night  the  second  order  mode  may  be  dominant  at  the  ground  to  distances  of 
several  thousand  kilometers.  For  the  first  mode  the  effect  of  height-gain  is 
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to  cause  the  dominance  of  the  second  node  to  be  greatly  diminished  or  even 
reversed*  This  height-gain  increase  with  height  for  the  first  order  node  has 
been  called  the  whispering  gallery  effect  by  Budden  and  Martin  (1962)  or  the 
earth-detached  node  by  Wait  and  Spies  (1963).  It  involves  successive  reflec¬ 
tions  from  the  ionosphere  with  essentially  no  intermediate  bounce  at  the 
ground. 

The  modal  constants  characterized  above  become  considerably  more  compli¬ 
cated  with  full  consideration  of  the  geomagnetic  field.  From  the  earliest 
days  of  VLF,  experimental  investigation  provided  evidence  of  an  apparent 
violation  of  the  reciprocity  principle.  Round  et  al  (1925),  for  example, 
found  that  VLF  radio  waves  apparently  suffered  significantly  larger  attenua¬ 
tion  for  propagation  in  a  generally  easterly  direction  them  in  the  reverse 
direction.  Some  very  early  results  presented  by  Wait  (1961)  for  sharply 
bounded  model  ionospheres  suggested  that  simple  harmonic  functions  may  ade¬ 
quately  describe  the  azimuthal  dependence  of  mode  constants  under  daytime 
conditions.  Further,  a  numerical  modeling  study  by  Ferguson  (1968)  showed 
that  for  daytime  conditions,  only  the  horizontal  conponent  of  the  geomagnetic 
field  which  is  transverse  to  the  direction  of  propagation  is  important.  This 
component  of  the  magnetic  field  is  the  source  of  nonreciprocity,  and  the 
remaining  components  were  commonly  assumed  only  to  alter  the  propagation 
characteristics  somewhat  (Makorov  et  al,  1970) .  The  most  severe  geomagnetic 
field  effects  were  thus  expected  in  the  equatorial  region  for  propagation 
across  the  magnetic  meridians. 

Perhaps  the  first  indications  of  extreme  complexity  in  modal  behavior  due 
to  the  geomagnetic  field  are  found  in  the  numerical  modeling  results  of  Snyder 
(1968a,  b)  and  Pappert  (1968).  Snyder  (1968a)  addressed  the  problem  of  mode 
numbering  for  typical  nighttime  ionospheres  as  a  function  of  azimuthal  varia¬ 
tion  for  a  horizontal  magnetic  field.  Results  indicated  that  a  mode  with 
lowest  west-to-east  phase  velocity  evolved  as  a  result  of  continuous  variation 
in  magnetic  azimuth  into  a  mode  in  the  reverse  direction  with  only  second 
lowest  phase  velocity.  Thus,  numbering  modes  in  terms  of  increasing  phase 
velocity,  as  is  commonly  done,  results  in  the  evolution  of  the  first  order 
west-to-east  mode  into  the  second  order  mode  in  the  opposite  direction,  and 
vice  versa.  Further  analysis  (Snyder,  1968b)  indicated  that  for  a  certain 
inclination  of  a  dipping  geomagnetic  field  and  at  a  specific  (near  north- 
south)  azimuth,  the  two  modes  considered  above  became  degenerate,  forming  a 
single  mode.  For  inclinations  more  nearly  horizontal  than  the  degenerate 
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conditions,  mode  numbering  inconsistencies  resulted;  but  for  more  nearly 
vertical  inclinations,  inconsistencies  in  mode  numbering  vanished.  Pappert 
(1968)  presented  the  results  of  numerical  modeling  for  easterly  propagation  at 
midlatitudes.  Assuming  a  vertical  dipole  exciter,  Pappert  found  that,  at  the 
high  frequency  end  of  the  VLF  band  (only  near  30  kHz),  inodes  which  are  princi¬ 
pally  transverse  electric  (TE)  may  be  of  importance  in  a  mode  sum.  Note  that 
a  TE  mode  can  be  excited  by  a  vertical  dipole  only  because  of  the  presence  of 
the  geomagnetic  field. 

Snyder  and  Pappert  (1969)  extended  the  analysis  to  consider  both  easterly 
and  westerly  midlatitude  nighttime  propagation  throughout  the  VLF  band  as  well 
as  to  include  azimuthal  dependencies  of  mode  parameters  for  a  central  VLF 
frequency  at  both  middle  and  equatorial  latitudes.  It  was  found  that  the 
importance  of  principally  TE  modes  is  much  more  pronounced  for  westerly 
propagation  at  midlatitudes  than  for  easterly  and  that  the  influence  of  the 
principally  TE  modes  on  the  mode  sum  can  be  significant  at  frequencies  as  low 
as  20  kHz  or  less  for  westerly  propagation.  Azimuthal  anomalies  included 
drastic  polarization  changes  in  going  from  easterly  to  westerly  paths.  For 
transverse  propagation  at  the  magnetic  equator,  it  was  shown  that  modes  which 
are  pure  transverse  magnetic  (TM)  for  propagation  to  the  east  may  be  pure  TE 
for  propagation  to  the  west.  This  is  tantamount  to  the  statement  that  modes 
which  have  dominant  excitation  for  easterly  propagation  may  have  vanishingly 
small  excitation  for  westerly  propagation.  Azimuthal  dependencies  are  charac¬ 
terized  by  rapid  variation  of  the  mode  constants  in  the  neighborhood  of  north- 
south  or  south-north  propagation,  for  both  equatorial  and  middle  latitudes. 
These  variations  manifest  themselves  in  marked  dif fei  ences  in  mode  sums  for 
azimuthal  changes  as  small  as  10°  or  less.  Further,  it  was  shown  that  the 
maximum  total  signal  attenuation  occurs  for  north-south  propagation  rather 
than  for  east-west  propagation. 

Although  Snyder  and  Pappert  (1969)  used  hypothetical  ionospheres,  their 
results  have  been  further  substantiated  by  Foley  et  al  (1973),  who  used 
"realistic  models  of  the  ionosphere."  Foley  et  al  essentially  duplicated  the 
analysis  of  Synder  and  Pappert  but  used  ionospheric  models  presented  in  the 
literature  by  Deeks  (1966)  and  Smith  et  al  (1968).  Foley  et  al  concluded  that 
mode  characteristics  are  dependent  on  the  ionosphere  model  assumed.  For 
example,  results  with  the  Smith  profile  were  similar  to  the  results  of  Snyder 
and  pappert,  whereas  under  certain  circumstances  results  differed  with  the 
Deeks  profile.  For  propagation  in  the  north-south  direction  with  the  Deeks 
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profile,  anomalies  were  found  in  the  behavior  of  modes  1  and  2.  It  would 
appear  that  between  16  and  20  kHz,  modes  1  and  2  interchange  their  identi¬ 
ties.  Although  unknown  to  Foley  et  al,  this  identity  interchange  would  seem 
to  be  just  another  manifestation  of  the  modal  degeneration  reported  by  Snyder 
( 1968b) . 
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THE  TERMINATOR 


A  manifestation  of  oultimode  propagation  is  found  in  the  phenomenon  of 
sunrise  and  sunset  fading,  which  is  characterized  by  periodic  and  repeatable 
variations  in  aiqplitude  and  phase  as  the  dawn-dusk  terminator  moves  along  a 
VLF  propagation  path.  Crombie  (1964)  explained  the  phenomenon  using  waveguide 
mode  concepts.  Of  all  possible  combinations  of  conditions  sunrise  and 
sunset,  position  of  terminator  on  propagation  path,  relative  positions  of 
transmitter  and  receiver,  etc  only  two  fundamental  situations  occur  within 
the  confines  of  Crombie* s  explanation. 

The  first  situation  depicts  the  field  incident  on  the  terminator  from  the 
transmitter  as  determined  by  two  modes,  and  the  field  at  the  receiver  as 
determined  by  only  one  mode.  This  situation  would  be  characteristic  of  a 
transmitter  on  the  night  side  of  the  terminator  and  the  receiver  on  the  day 
side  sufficiently  removed  to  avoid  multimode  propagation.  The  received  signal 
then  depends  on  the  amplitude  and  phase  relationships  of  the  modes  incident  on 
the  terminator  and  on  the  efficiency  with  which  these  modes  are  converted  by 
the  terminator  to  the  single  mode  which  reaches  the  receiver.  Because  the 
modal  relationships  are  determined  by  the  distance  of  the  terminator  from  the 
transmitter,  periodic  variations  will  occur  in  the  received  signal  as  the 
terminator  moves  along  the  propagation  path.  The  periods  will  be  determined 
by  the  difference  in  the  phase  velocities  of  the  two  modes  entering  the 
terminator,  and  the  variations  must  occur  simultaneously  for  all  points  in  the 
portion  of  the  propagation  path  beyond  the  terminator  region.  The  fading 
would  be  most  intense  when  the  two  entering  modes  have  the  most  pronounced 
interference  ~  usually  when  the  terminator  is  near  the  transmitter. 

The  second  situation  depicts  the  incident  field  as  a  single  mode.  The 
presence  of  the  terminator  leads  to  multimode  generation  in  the  region  past 
the  terminator,  which  in  turn  leads  to  modal  interference.  In  this  situation, 
the  interference  pattern  would  appear  to  be  attached  to  the  terminator  and  to 
move  synchronously  with  it.  Signal  minima  are  then  to  be  ejected  at  fixed 
distances  from  the  terminator,  the  distance  between  fades  being  determined  by 
the  difference  between  the  phase  velocities  o:  the  two  modes  leaving  the 
terminator  region.  The  deepest  signal  fades  would  be  expected  as  the  termin¬ 
ator  approaches  the  receiver.  This  situation  corresponds  to  propagation  from 
the  sunlit  side  of  the  terminator  into  the  night  side.  Note  that  in  both 
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situations,  the  fade  spacing  would  be  determined  by  the  differences  between 
the  phase  velocities  of  two  nighttime  nodes* 

A  thorough  experimental  examination  of  the  Crombie  model  was  conducted  by 
Walker  (1965).  who  recorded  16  kHz  transmissions  from  NBA  Balboa,  Panama, 
while  making  repeated  shipboard  crossings  of  the  Atlantic  in  the  equatorial 
region.  During  sunrise  (propagation  from  night  to  day)  signal  fades  were 
observed  simultaneously  at  all  points  in  the  day  portion  of  the  path  in 
agreement  with  the  first  situation  discussed  above.  Walker  determined  the 
change  in  terminator  position  between  successive  fades  and  from  this  calcu¬ 
lated  the  required  differences  in  the  phase  velocities  of  the  two  interfering 
modes.  These  results  were  in  good  agreement  with  theoretical  results  for 
nighttime  ionospheres,  further  substantiating  the  Crombie  model.  Observations 
were  also  made  during  sunset  (propagation  from  day  to  night)  and  results  were 
in  excellent  agreement  with  the  Crombie  model. 
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EQUATORIAL  ANOMALY 


The  observations  just  discussed  were,  for  the  most  part,  made  at  middle 
or  low  latitudes.  As  the  previous  discussion  indicates,  such  observations 
have  been  explained  adequately  by  the  use  of  quite  single  VLF  waveguide  mode 
theory.  However,  some  observations  of  VLF  transmissions  over  long  transequa- 
torial  paths  have  apparently  indicated  an  anomalous  effect  associated  with 
nighttime  or  transitional  transequator ial  VLF  propagation. 

Of  the  reports  on  transequatorial  VLF  propagation,  one  of  the  earliest  to 
show  an  apparent  anomaly  was  by  Chilton  et  al  (1964).  They  made  simultaneous 
observations  of  NBA  transmissions  (Canal  Zone,  Panama)  at  18  kHz  for  paths  to 
both  the  northern  and  southern  hemispheres.  Both  paths  were  of  similar 
length,  but  the  southern  hemisphere  path  crossed  the  magnetic  equator. 
Although  they  did  not  look  at  transition  fading,  they  did  observe  an  anomalous 
difference  in  the  diurnal  changes  of  signal  amplitude  and  phase.  Chilton  et 
al  suggested  that  their  observations  resulted  from  a  difference  in  ionization 
profile  due  to  latitudinal  dependence  of  cosmic  rays.  More  recently,  Chilton 
and  Crary  (1971)  suggested  that  the  latitude-dependent  ionization  source  might 
arise  from  the  recently  discovered  x-ray  stars. 

Araki  et  al  (1969)  recorded  dual  frequency  (15.5  and  22.3  kHz)  VLF  sig¬ 
nals  in  Japan  from  station  NWC  in  Australia.  They  reported  sunrise  fading  on 

22.3  kHz  but  not  on  15.5  kHz  and  further  noticed  that  during  fading  at  22.3 
kHz  the  first  phase  change  was  an  increase  followed  by  a  phase  decrease.  The 
nighttime  field  strength  was  somewhat  less  at  15.5  kHz  and  somewhat  higher  at 

22.3  kHz  than  during  the  day.  Araki  et  al  (1969)  concluded  that  their  results 
suggest  the  existence  of  a  frequency-dependent  anomaly  of  the  nighttime  field 
relative  to  the  daytime  field  along  their  transequatorial  propagation  path. 

Bickel  et  al  (1970)  conducted  an  experiment  to  investigate  the  apparently 
anomalous  results  of  Snyder  and  Pappert  (1969),  who  predicted  maximum  signal 
attenuation  in  a  southerly  rather  than  westerly  direction.  Measurements  of 

23.4  kHz  signals  from  station  NPM  (Lualualei,  Hawaii)  were  made  aboard  an 
airplane  as  it  flew  radials  along  the  propagation  paths  to  Seattle,  Ontario 
(California),  Samoa,  and  Wake  Island.  Bickel  et  al  confirmed  the  Snyder  and 
Pappert  (1969)  results  and  concluded  that  the  increased  southerly  attenuation 
is  an  effect  of  the  geomagnetic  field.  An  unexplained  effect,  however,  was 
observed.  For  the  radial  to  Samoa,  rapid  signal  amplitude  variations  with 
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distance  were  observed  in  the  vicinity  of  the  geomagnetic  equator.  Although 
equipment  problems  may  be  suspected,  other  Hawaii-based  transmissions  simul¬ 
taneously  recorded  by  means  of  completely  separate  equipments  displayed  simi¬ 
lar  results  (Bickel,  private  communication,  1970). 

The  first  report  of  an  equatorial  anomaly  associated  with  transition 
fading  was  made  by  Lynn  (1967).  He  reported  on  VLF  transmissions  at  18.6  kHz 
from  NLK  (Jim  Creek,  Washington)  to  Smithfield,  South  Australia.  This  propa¬ 
gation  path  is  approximately  13  500  km  long,  and  the  receiver-to-geomagnetic- 
equator  distance  along  the  path  is  approximately  5700  km.  The  direction  of 
propagation  is  essentially  southwesterly.  The  interference  distance  reported 
by  Lynn  for  sunrise  transition  fading  is  approximately  2000  km  for  the  termin¬ 
ator  located  in  midlatitudes  (in  excess  of  ±20*  from  the  geomagnetic  equator). 
Whenever  sunrise  transition  fading  was  observed  while  the  terminator  was 
within  ±20*  of  the  geomagnetic  equator,  the  interference  distance  increased  to 
as  much  as  3700  km;  and  when  averaged  over  the  anomaly,  it  had  a  value  of  2900 
km.  Lynn  (1967)  concluded  from  his  observations  that  the  change  in  interfer¬ 
ence  distance  resulted  from  a  change  in  the  difference  of  phase  velocity  for 
two  modes  as  well  as  from  a  change  in  the  relative  phases  of  the  appropriate 
mode  conversion  coefficients.  He  left  unanswered  the  question  of  a  possible 
cause  for  these  changes.  It  is  interesting  to  note  that  Lynn's  conclusions 
are  dependent  on  his  selection  of  both  a  VLF  propagation  model  and  a  mode  con¬ 
version  model.  Lynn  (1967)  assumed  an  isotropic  propagation  model  wherein  any 
changes  in  propagation  parameters  necessitates  changes  in  ionospheric  condi¬ 
tions. 

Kaiser  (1968)  examined  VLF  signals  at  18.6  and  24  kHz  from  NLK  and  at  20 
kHz  from  WWVL,  received  at  Lower  Hutt,  New  Zealand.  These  propagation  paths 
are  similar  to  the  path  considered  by  Lynn  (1967),  and  Kaiser's  results  are 
much  the  same  as  Lynn's.  An  anomalously  large  interference  distance  (in 
excess  of  3000  km)  was  observed  for  sunrise  transition  fading  when  the  termin¬ 
ator  was  within  about  ±20*  of  the  magnetic  equator.  Kaiser  also  examined 
transition  fading  data  previously  published  elsewhere  and  concluded  that 
anomalously  large  values  of  the  sunrise  interference  distance  are  characteris¬ 
tic  of  generally  east-to-west  VLF  propagation  below  about  30*  magnetic  lati¬ 
tude.  A  further  conclusion  by  Kaiser  (1968)  is  that  any  anomaly  in  sunset 
transition  fading  for  such  paths  is  relatively  minor.  Kaiser  suggests  that  a 
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possible  cause  of  the  anomaly  might  be  a  larger  and  sharper  day-night  wave¬ 
guide  transition  near  the  geomagnetic  equator  as  compared  to  higher  latitudes. 
He  does  not  suggest,  however,  a  possible  cause  for  the  latitudinal  dependence 
of  the  day-night  waveguide  transition. 

Additional  observations  of  this  transition  fading  anomaly  have  been 
reported  by  Lynn  (1969,  1970).  These  observations  further  substantiate  the 
general  characteristics  of  the  anomaly.  A  theoretical  interpretation  of  the 
transequator ial  anomaly  was  presented  by  Lynn  (1970).  His  interpretation  is 
based  only  on  assumed  latitudinal  variations  in  modal  phase  velocities  and 
relative  phases  of  the  mode  conversion  coefficients.  He  does  not  suggest  a 
possible  cause  for  the  variations. 

Meara  (1973)  has  extended  Lynn's  (1970)  analysis  to  include  determination 
of  the  changes  in  phase  velocities  within  the  equatorial  anomaly  region.  He 
concludes  that  the  phase  velocity  of  the  first  mode  is  essentially  unaltered 
by  propagation  through  the  equatorial  region,  whereas  that  of  the  second  mode 
is  reduced.  Note  that  Meara' s  results  are  entirely  dependent  'on  the  assumed 
propagation  and  mode  conversion  models.  Once  again,  no  explanation  of  a 
possible  cause  for  the  variations  in  phase  velocities  is  given. 

A  further  application  of  Lynn's  (1970)  analysis  to  explain  some  anomalous 
diurnal  changes  in  transequatorial  VLF  propagation  has  been  made  by  Araki 
(1973).  The  analysis  follows  essentially  that  previously  discussed.  Araki 
does  go  a  step  further  than  Lynn  (1970)  and  Meara  (1973),  however,  in  that  he 
assumes  the  validity  of  the  isotropic  propagation  model  and  asserts  that  any 
change  in  propagation  parameters  results  from  a  modification  of  the  nighttime 
equatorial  ionosphere.  Lynn  (1978)  has  pointed  out  that  the  propagation  paths 
of  Chilton  et  al  (1964)  and  Araki  (1973)  have  a  west-east  exponent  rather 
than  an  east-west  component  as  have  the  other  paths  discussed  here.  He 
further  claims  that  only  a  minor  azimuthal  variation  in  relative  phase 
velocities  of  interfering  nighttime  modes  can  account  for  the  observations  of 
Chilton  et  al  (1964)  and  Araki  (1973).  Lynn  (1978)  thus  concludes  that  the 
observations  of  Chilton  and  Crary  (1971)  cannot  be  seen  as  direct  evidence  for 
the  control  of  night  VLF  reflection  heights  by  stellar  X-ray  sources.  But 
Svennesson  and  Westerlund  (1979),  using  waveguide  mode  theory  combined  with 
nighttime  profiles  given  by  ionospheric  theory,  claim  results  that  support  the 
reports  of  X-ray  star  effects  on  VLF  radio  signals. 


To  address  the  question  of  transequatorial  VLF  propagation  at  night  or 
Airing  transition  periods  requires  theoretical  VLF  radio  wave  propagation 
Modeling  that  includes  consideration  of  inhoiaogeneities  along  the  propagation 
path  as  well  as  in  the  vertical  direction.  Anisotropic  effects  due  to  the 
geomagnetic  field  mist  also  be  considered*  A  propagation  model  with  these 
capabilities  has  been  developed  and  is  discussed  in  the  following  sections. 
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2.  WAVEGUIDE  NODE  THEORY  -  BOKiaOMTALLY  HOHOGEHBOOS  GUIDES 


COORDINATE  TRANSFORMATION 

A  review  of  the  waveguide  mode  theory  for  VLF  propagation  in  the  earth- 
ionosphere  waveguide,  homogeneous  in  and  transverse  to  the  direction  of 
propagation,  is  presented  in  this  section.  The  earth-ionosphere  guide  is 
modeled  as  cylindrical,  with  effects  of  curvature  transverse  to  the  direction 
of  propagation  accounted  for  through  the  use  of  an  appropriate  "spreading 
term"  as  described  in  Section  1  Also,  as  mentioned  there,  the  formulation 
enployed  in  this  study  is  based  on  the  pioneering  work  of  Budden  (1961, 
1962).  Because  Budden's  formulation  is  essentially  planar,  a  coordinate 
transformation  is  employed  to  convert  between  cylindrical  geometry  and  planar 
geometry.  The  transformation,  introduced  by  Richter  (1966)  with  respect  to 
earth  flattening  techniques  and  used  by  Pappert  (1968)  and  Abbas  et  al  (1971) 
in  VLF  propagation  studies,  is  briefly  outlined. 

Application  of  the  conformal  coordinate  transformation  of  the  form  p  * 
aez/a;  $  ■  x/a  results  in  the  conversion  of  circles  (p  *  constant)  in  the  p,$ 
plane  to  straight  lines  (2  *  constant)  in  the  x,z  plane.  If  a^ represents  the 
radius  of  the  earth,  the  transformation  maps  the  surface  of  the  earth  (P  =  a) 
to  the  line  2=0.  The  transformation  further  maps  radii  ( 4>  *  constant)  from 
the  p,<J>  plane  to  lines  x  =  constant  in  the  x,z  plane.  A  full  circular  arc  (♦ 
=  2ff)  at  P  =  0  is  mapped  so  that  the  length  of  x  is  the  circumference  of  the 
earth. 

The  transformation  used  here  differs  from  that  enployed  by  Richter 
(1966),  Pappert  (1968),  and  Abbas  et  al  (1971)  in  that  they  maintained  some 
selected  altitude  (h)  as  invariant  in  both  coordinate  frames.  However,  the 
two  transformations  differ  only  in  quantities  of  second  order  in  nallness. 
As  the  invariant  altitude  (h)  approaches  zero,  the  two  transformations  become 
identical. 

Maxwell's  curl  equations  with  time  variation  of  the  form  elwt  for  a 
medium  characterized  by  a  dielectric  tensor  may  be  written  as 


V  x  £  -  -ikoa 

l  *  «  “  ikDS  *  V 


(3) 
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where  kQ  is  the  wave  number  of  free  apace  (k0  »  w/c)  and  £  is  the  dielectric 
tensor.  In  equations  (3),  the  variable  JJ  is  actually  the  magnetic  intensity 
multiplied  by  the  free  space  iopedance.  Assuming  3/3y  *  0,  equation  (3)  in 
cylindrical  coordinates  may  be  written 
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Applying  the  coordinate  transformation  to  equation  (4): 
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By  making  the  identifications 


Ex  =  p/a  fy  Ey  =  V  Ez  =  P/a  V 

Hx  =  p/a  H^i  Hy  =  Hy}  Hz  =  p/a  Hp,  (6) 


equations  (5)  become  identical  to  what  would  be  obtained  in  a  rectangular 
coordinate  frame  for  a  medium  characterized  by  permeability  tensor  and 
dielectric  tensor  of  the  form 
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Some  additional  insight  is  obtained  by  considering  the  transformed 
equations  in  the  free  space  region  between  the  earth  and  the  ionosphere. 
Making  the  usual  substitutions  into  equation  (5c)  from  equations  (5a)  and 
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(5b),  along  with  the  free  space  assumption,  yields  the  equations 
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(9a) 
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Thus,  the  effect  of  the  transformation  is  to  replace  the  original  cylindrical 
free  space  region  with. a  planar  region  filled  with  a  medium  which  has  refrac¬ 
tive  index  equal  to  eZ//a  for  TE  polarized  waves  or  which  has  permeability 
equal  to  for  TM  polarized  waves, 

MODE  EQUATION 

The  prototype  waveguide  model  considered  here  consists  of  two  plates  of 
infinite  extent,  parallel  to  the  x-y  plane  and  located  at  z  »  0  and  z  -  h. 
The  interior  of  the  guide  is  assumed  to  be  free  space.  This  prototype  wave¬ 
guide  corresponds  to  the  "fictitious  free  space  region"  to  be  later  discussed 
in  this  section,  where  it  is  used  to  determine  modal  excitation  factors  and 
polarization  coupling.  The  electromagnetic  field  is  assumed  to  be  independent 
of  the  y-coordinate  and  propagation  is  assumed  to  be  in  the  x-direction.  Any 
given  waveguide  mode  is  considered  to  be  composed  of  two  crossing  plane  waves, 
upgoing  and  downgoing,  with  wave  normals  at  the  angles  ±0  to  the  x-axis. 

The  nature  of  the  boundary  planes  will  be  specified  in  terms  of  reflec¬ 
tion  coefficients.  These  coefficients  are  defined  in  terms  of  the  ratios  of 
specific  components  of  the  upgoing  and  downgoing  waves  and  are  therefore 
dependent  on  the  polarization  of  the  waves.  A  convenient  way  of  describing 
the  polarization  is  in  terms  of  the  linear  components  of  the  complex  electric 
vector  which  are  perpendicular  and  parallel  to  the  plane  of  incidence  (the  x-z 
plane).  The  notation  1(1)  will  indicate  that  the  electric  vector  is  parallel 
(perpendicular)  to  the  plane  of  incidence. 
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In  general,  the  upper  boundary  nust  be  considered  anisotropic.  Thus 
there  can  be  a  change  in  polarization  upon  reflection.  The  reflecting 
properties  of  the  upper  boundary  are  given  by  reflection  coefficients  having 
the  matrix  form 


B  = 


(10) 


where  the  first  subscript  describes  the  polarization  of  the  incident  (upgoing) 
wave  and  the  second  subscript  denotes  the  polarization  of  the  reflected  wave. 
The  properties  of  the  lower  boundary  are  assumed  isotropic  and  specified  by  a 
diagonal  reflection  matrix  of  the  form 


R 

M 


(11) 


If  the  upgoing  wave  is  considered  as  the  primary  or  incident  wave,  the 
downgoing  wave  can  be  considered  as  resulting  from  reflection  at  the  upper 
boundary.  This  wave  will  in  turn  be  reflected  at  the  lower  boundary,  giving 
rise  to  another  upgoing  wave.  The  requirement  that  the  second  upgoing  wave  be 
the  same  as  the  original  yields  the  mode  condition.  This  process  is  shown  in 
the  diagram  below. 


The  primary  wave  is  represented  Ify  the  column  matrix  u  .  The  second  upgoing 

•o 

wave  is  therefore  given  by  the  relationship 


u  «  R  R  u  exp[-2ikhsin( 6)) 


where  the  exponential  term  accounts  Cor  phase  change  in  the  double  passage 
across  the  guide.  The  self-consistency  requirement  for  the  waveguide  mode  is 
thus 


det  {j  g  exp[-2ikhsin( 9 ) ]  -  =*  0, 


(12) 


where  1  is  the  unit  diagonal  matrix. 


POINT  SOURCES  AND  LINE  SOURCES 

The  source  is  assumed  to  be  a  vertical  Hertzian  dipole.  Budden  (1961 , 
1962)  has  shown  that  a  vertical  dipole  is  equivalent  to  an  infinite  set  of 
line  quadrupoles,  where  the  line  quadrupoles  are  located  parallel  to  the  x-y 
plane  and  inclined  at  (congjlex)  angles  to  the  y-axis.  This  is  easily  seen  by 
considering  the  Hertz  vector  for  a  dipole  at  the  origin.  If  the  dipole  is 
parallel  to  the  z-axis,  the  Hertz  vector  has  only  a  z-component  and  is 
proportional  to  exp(-ikr)/r .  For  a  dipole  of  unit  strength,  the  Hertz  vector 
is  as  follows  (Budden,  1962): 
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-ikr 
4  nr 


k_  (■  if/2+i« 
8*  J-ir/2-i“ 


cos(a) 


(ks)  dot. 


(13) 


where 


a  *  arctan  {z/[x  cos(w)  +  y  sin(<u)]} 

>j2  2  1 1/2 

s  *  {[x  cos(ttf)  +  y  sin{w) J  +  *  ) 


Thus  a  is  the  distance  from  a  point  (x,  y,  z)  to  a  line  in  the  x-y  plane 
passing  throu^i  the  origin  at  an  angle  u  to  the  y-axis,  and  a  is  the  angle 
between  the  x-y  plane  and  the  line  along  s.  is  a  Hankel  function  of  the 
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second  kind,  of  order  one.  The  term  cos(  a)H^  ( 2 )  ( ^g)  would  result  from  a  line 
quadrupole  source  in  the  x-y  plane  at  an  angle  w  to  the  y-axis.  Consider  a 
lint'  quadrupole  source  as  composed  of  two  line  dipole  sources  at  x  *  and 

x  =  -  V2  «X  in  the  x-y  plane.  The  line  dipoles  consist  of  closely  spaced 
dipoles  oriented  parallel  or  antiparallel  to  the  z-axis  and  distributed  paral¬ 
lel  to  the  y-axis.  The  Hertz  vector  for  the  line  dipoles  is  proportional  to 
H  <2)(kp)  where  p  =  (x2  +  z2)1//2.  Letting  6x  tend  to  zero  and  letting  the 
strengths  (M,  -M)  of  the  line  dipoles  increase  in  such  a  way  that  the  product 
Mfix  remains  finite  will  result  in  a  term  proportional  to  3hq  (kp)/3x.  which 
reduces  to  a  term  of  the  form  cost  9)  2^  (kp ) ,  where  cos(9)  =  x/p.  From  the 

definitions  of  s  and  a  in  (14),  it  is  clear  that  the  integral  in  (13)  results 
from  an  infinite  set  of  line  quadrupoles  in  the  x-y  plane,  each  at  a  different 
angle  u)  to  the  y-axis.  Because  the  path  of  integration  in  (13)  is  not  unique, 
there  is  an  infinite  number  of  ways  in  which  the  angles  w  may  be  chosen. 

Consider  a  guide  bounded  by  perfectly  conducting  plates  with  a  line 

source  parallel  to  the  y-axis.  This  source  generates  waveguide  modes,  and 

each  consists  of  crossing  plane  waves  with  normals  in  the  x-y  plane  at  angles 

±6  to  the  x-axis.  For  a  line  source  at  an  angle  w  to  the  y-axis,  the  com- 
m 

ponent  plane  wave  has  the  form  exp(-ik{[x  cos(w)  +  y  sin(u>)]Cm  +  z  SjJJ*  where 

C  «  cos(  0  )  and  S  =  sin(  0_) .  A  point  source  within  the  guide  is  now  repre- 
m  mm  m 

sented  by  an  infinite  set  of  line  sources  such  that  the  total  field  of  the 
component  wave  is  proportional  to 


/  ^2+i»  exp(-ik{  [x  cos(w)  +  y  sin(u)]cm  t  z  S^} )  d<u  .  05) 


Making  the  substitutions  p 


(x2  +  y2)1^2  and  cos(u)  =  x/P  in  (15), 


/-x/2-i»  exp(_lkpcmcos(  “~u> )  ejq>(-ikzSa)  do>  . 


06) 
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The  integration  in  (16)  yields  a  Hank el  function  of  the  second  kind,  of  order 
zero  (Sommerfeld,  1949).  Thus  the  total  field  of  the  coiqponent  wave  is 
proportional  to 


Ho<2>(kpc»)  exp(-ikzsm). 

The  appropriate  Hankel  function  exchange  for  the  exponential  term  when  a  line 
source  is  used  is  therefore 


exp(-ikxcj  -  H0<2>(kPCn). 


(17) 


MODAL  EXCITATION  FACTORS  AND  POLARIZATION  COUPLING 

The  excitation  of  waveguide  inodes  is  conveniently  treated  in  terms  of 
line  sources.  The  technique  consists  of  selecting  one  of  the  possible  line 
sources  for  the  integrand  in  (13),  finding  the  mode  excitation  for  this  line 
source,  then  effecting  a  change  from  an  exponential  term  to  the  appropriate 
Hankel  function.  Following  the  method  of  Budden  (1961,  1962),  the  source  is  a 
line  of  quadrupoles  at  height  d  (near  the  ground)  parallel  to  the  y-axis.  A 
fictitious  free  space  region  is  introduced  that  extends  from  z  =  1  (lower)  to 
z  -  u  (upper)  such  that  1  <  d  <  u  (the  line  source  is  assumed  within  the  free 
space  region) .  The  region  above  z  =  u  is  characterized  by  the  reflection 
coefficient  matrix  R^,  while  the  region  below  z  *  1  is  characterized  by  a 
reflection  coefficient  matrix  §^.  If  the  line  source  has  unit  strength,  the 
yconponent  of  the  magnetic  field  in  infinite  free  space  can  be  written 

k  3 

■y-  4if"  /0:;r  exp(iko(xCHz-d|S))c2d8  ,  (18) 

1  o 

where  C  *  cost 8),  s  -  sin(0),  and  the  limits  of  integration  are  selected  to 
ensure  convergence  for  all  x  ( Somnerfeld,  1949).  The  primary  field  integrand 
in  equation  (18)  represents  a  plane  wave  which  undergoes  reflections  at  z  >  u 
and  z  -  1,  and  these  reflected  waves  must  be  added  to  the  primary  field  to 
obtain  the  total  field.  Four  combinations  are  possible,  two  from  waves 
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initially  upgoing  or  downgoing  at  the  source  and  two  from  waves  reaching  the 
reception  point  upgoing  or  downgoing. 

\n  initially  upgoing  wave  is  converted  to  a  downgoing  wave  by  the  reflec¬ 
tion  process  denoted  by  the  operator  (g  n  and  to  an  upgoing  wave  by  the 

operator  (g  g  J  .  An  initially  downgoing  wave  is  converted  to  an  upgoing 

**■  ^  . _  -  M 

wave  by  the  operator  g^  and  to  a  downgoing  wave  by  the  operator 

(r  r  )n  .  The  N  denotes  N  reflection  processes  of  the  the  type  described  by 
"U"l 

the  bracketed  term.  With  each  reflection  process  mist  also  be  associated 
phase  terms  to  account  for  the  transversals  to  and  fro  between  z  »  1  and  z  « 
u.  The  subscripts  u  and  1  denote  reflection  coefficients  evaluated  at  z  »  u 

and  z  =  1,  respectively. 

The  total  field  is  obtained  by  combining  the  primary  and  all  reflected 
fields  and  letting  the  dimensions  of  the  free  space  region  tend  to  zero.  This 
requires  taking  the  limits  1  ♦  u  and  z  +  d  in  such  a  way  that  1  <  d  <  u  is 
maintained.  The  total  field  may  be  written  as 
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The  reflection  coefficients  are  evaluated  at  i  *  d  and  the  subscripts  1  and  u 
are  no  longer  needed.  in  equation  (19),  the  first  summation  starts  with  N»0 
in  order  to  include  the  primary,  nonreflected  wave;  and  g  is  a  column  matrix 
whose  second  element  is  zero,  corresponding  to  a  vertical  electric  dipole 
source. 

The  order  of  summation  and  integration  is  now  reversed  so  that  each 
integrand  will  contain  a  geometric  series.  It  can  be  shown  (Budden,  1961, 
1962)  that  the  series  converges  for  certain  simplified  cases,  and  convergence 
is  assumed  here.  Summation  of  the  geometrical  series  gives 
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Consequently  the  result  of  summing  the  series  in  equation  (19)  is 
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If  the  sign  ,of  0  is  reversed,  the  role  of  incident  and  reflected  wave  is 
exchanged  and  each  reflection  element  is  replaced  by  its  inverse.  The 

integrand  simply  changes  sign  upon  reversal  of  the  sign  of  6.  The  integral 

may  be  evaluated  by  using  the  theory  of  residues,  deforming  the  contour  to  run 

symmetrically  through  the  origin  from  it/2  +  i°»  to  -*/2  -  i°°.  Because  the 
integrand  is  antisymmetric  about  0=0,  the  integral  along  the  new  contour  is 
zero.  In  changing  the  contour,  however,  singularities  are  crossed.  The  most 
important  of  these  are  poles  of  the  factor  (l  -  £  g)_1  •  The  residue 

contribution  of  a  simple  pole  at  @n  is  as  follows: 
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where 


F  -  det  (l  -  «  ft) 
x  »  lim  [fCx  -  5  S)"1], 

n  0-*-®n 


and  the  subscript,  n,  signifies  that  the  quantities  are  evaluated  at  0  =  0ft. 
The  pole  condition  F  =  0  is  identical  to  the  waveguide  mode  equation  (12)  in 
the  limit  h  ♦  0.  Thus  the  solution  of  the  pole  condit  .on  yields  the  waveguide 
modes,  and  the  residues  at  the  poles  give  the  excitation  factors.  Performing 
the  matrix  multiplications  in  equation  (24),  using  the  Hankel  function  substi¬ 
tution  [equation  (17)],  and  employing  the  geometrical  spreading  factor  yields 
the  equation 
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To  be  compatible  with  the  mode  sun  of  Section  1  [equation  (1)],  it  is  neces¬ 
sary  to  change  the  reference  direction  for  0  from  the  x-axis  to  the  z-axis. 
Doing  so  simply  requires  interchanging  the  factors  S  and  sin(9)  with  the  fac¬ 
tors  C  and  cos( 6 ) . 

Equation  (26)  requires  the  source  and  reception  point  to  be  at  the  same 
height,  which  is  also  the  height  at  which  the  reflection  elements  are  to  be 
determined.  If  the  dependence  of  on  height  z  is  given  by  the  function 
Gn(z),  then  an  obvious  extension  of  equation  (26)  to  include  arbitrary  source 
and  reception  height  involves  multiplication  by  the  function 


Gn(ZT)  Gn^R^ 


G  (d)  G  (d) 
n  n 


> 


where  is  the  transmitter  height,  Zjj  is  the  receiver  height,  and  d  is  the 
evaluation  height  for  the  reflection  elements,  assumed  to  be  below  the 
ionosphere. 

Whenever  the  ionosphere  is  assumed  to  be  isotropic  (no  geomagnetic 
field),  the  field  components  can  be  divided  into  two  sets.  These  are  the 
well-known  TE  (Hy  -  0)  and  TM  (Ey  =  0)  modes.  When  anisotropy  is  included, 
however,  a  given  mode,  in  general,  involves  both  Hy  and  Ey.  To  determine  the 
coupling  between  the  TM  and  TE  components  of  the  eigenfunctions,  we  again 
introduce  a  thin  fictitious  region  of  free  space  within  the  transformed  planar 
waveguide.  The  region  is  assumed  to  be  centered  at  some  altitude  z  «*  d  which 
is  below  the  ionosphere,  so  that  the  region  z  <  d  can  be  considered  isotropic. 
Within  this  region,  the  fields  of  the  TM  and  TE  conponents  can  be  considered 
in  terms  of  upgoing  and  downgoing  plane  waves.  Furthermore,  the  eigenfunc¬ 
tions  satisfy  the  condition  that  a  particular  component,  say  the  upgoing,  must 
be  returned  to  its  original  value  after  being  reflected  from  the  two  bound¬ 
aries  of  the  thin  free  space  region.  Thus,  for  the  upgoing  component  e“  (the 
component  perpendicular  to  the  plane  of  propagation) ,  the  requirement  is  that 
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or 


In  the  free  space  region 
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cos(6)  eu  -  sin(9)  eu  =  hU 
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where  is  the  upgoing  component  of  Hy  and  e^  is  the  upgoing  component  of 
Ey.  The  total  Hy  and  Ey  fields  within  the  free  space  gap  at  z  =  d  are 
therefore  represented  by  the  equations 
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and  their  ratio  is  as  follows: 
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In  equation  (28)  the  reflection  coefficients  ere  evaluated  at  *  «  d. 
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If  significant  variations  occur  (as  they  often  do)  along  the  direction  of 
propagation  in  the  earth- ionosphere  waveguide,  then  homogeneous  guide  formula¬ 
tions  require  modification.  For  such  cases,  there  are  two  techniques  for  mode 
summation  that  are  commonly  used  at  VLF.  One  employs  a  WKB  approach;  the  other 
uses  mode  conversion. 

In  the  WKB  technique,  an  eigenmode  is  assumed  to  be  uniquely  and  indepen¬ 
dently  identifiable  anywhere  along  a  propagation  path.  Furthermore,  each  mode 
is  assumed  to  depend  only  on  the  local  characteristics  of  the  waveguide  and  to 
propagate  independently  of  the  existence  of  any  other  modes.  The  WKB  formula¬ 
tion  generalizes  the  expression  for  the  mode  summation  in  equation  (1)  to  the 
form 
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Note  that  the  only  differences  arise  from  (1)  the  use  of  the  geometric  mean  of 
the  excitation  factors  at  transmitter  (AT)  and  from  receiver  ( AR]  locations 
and  (2)  the  integration  of  the  complex  propagation  factor  along  the  path.  The 
WKB  formulation  was  apparently  introduced  simply  as  an  "obvious  extension"  of 
the  homogeneous  mode  sum  by  Wait  (1964),  who  further  noted  that  the  form  is 
similar  to  what  would  be  expected  if  a  WKB  formulation  were  applied,  although 
no  actual  derivation  was  presented.  The  use  of  the  name  WKB  mode  sum  has 
persisted. 

In  the  mode  conversion  method,  eigenmode s  are  assumed  to  be  uniquely  and 
independently  identified  only  in  a  local  sense.  Even  though  the  eigen  charac¬ 
teristics  of  a  mode  depend  only  on  local  characteristics  of  the  waveguide,  the 
propagation  of  each  mode  depends  not  only  on  the  existence  of  other  modes  but 
also  on  the  past  and  future  propagation  history  of  itself  and  the  other  modes. 
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The  approach  that  uses  mode  conversion  does  not  lend  itself  conveniently 
to  consideration  of  transition  fading  on  very  long  VLF  propagation  paths.  The 
WKB  method,  although  easily  applied,  is  unable  to  handle  the  rapid  change  in 
guide  characteristics  that  occurs  at  the  terminator  (Pappert  and  Snyder, 
1972) .  A  difficulty  encountered  when  using  mode  conversion  arises  from  the 
necessity  of  knowing  the  future  propagation  history  of  each  mode.  Making 
appropriate  assumptions  to  relieve  this  requirement  can  lead  to  some  simpli¬ 
fication,  but  the  resulting  formulation  is  still  not  convenient  for  very  long 
paths.  A  further  difficulty  is  that  conventional  mode  conversion  techniques 
require  a  knowledge  of  all  the  vector  components  of  the  electromagnetic  field 
throughout  all  space. 

Numerous  methods  have  been  presented  for  circumventing  the  latter 
problem.  Halt  ( 1968a,  b) ,  for  exasple,  formulated  mode  conversion  by  using 
isotropic  impedance  boundary  conditions,  with  the  result  that  fields  were 
required  only  for  the  region  between  the  upper  (ionospheric)  and  the  lower 
(ground)  impedance  boundaries.  Such  boundary  conditions  are  certainly  not 
appropriate  for  the  nighttime  ionosphere  and  thus  for  consideration  of  trans¬ 
itions  from  night  to  day.  Smith  (1974)  extended  Wait's  (1968a,  b)  results  to 
include  boundary  anisotropy  by  using  inpedance  matrix  boundary  conditions. 
His  analysis  however,  was  for  a  planar,  empty  waveguide  and  was  directed 
toward  considering  ground  effects. 

Galejs  (1971)  presented  a  modified  mode  conversion  formulation  for  an 
abrupt  change  between  daytime  and  nighttime  anisotropic  ionospheres.  Conven¬ 
tional  mode  conversion  theory  usually  represents  the  EM  field  as  an  expansion 
in  terms  of  a  suitably  chosen  complete  set  of  orthogonal  functions.  The 
Galejs  formulation  enploys  an  inconplete  set  of  nonorthogonal  functions.  The 
formulation  now  to  be  employed  is  similar  to  that  of  Galejs  (1971). 

OPERATOR  EQUATIONS 

The  operator  equations  presented  in  this  section  are  based  on  Maxwell's 
equations  in  two  (rectangular)  dimensions  and  are  patterned  after  the  formu¬ 
lations  of  Pappert  and  Smith  (1972).  Furthermore,  material  media  are  assumed 
to  be  characterised  by  constitutive  parameters  (permittivity  and  permeability) 
that  can  be  inhomogeneous  in  one  dimension.  This  is  quite  compatible  with  the 
usual  procedure  for  handling  the  earth- ionosphere  waveguide  by  describing  the 
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electromagnetic  properties  <>f  the  ionosphere  with  the  aid  of  a  vertically 
inhomogeneous  dielectric  tensor. 

We  start  with  equation  (5),  use  the  substitutions  in  equation  (6),  and 
eliminate  the  Ex  and  Hx  components.  This  leads  to  the  differential  operator 
equation 


^  ~  ik  3x‘ 
o 

where  e  is  a  four-element  column  vector: 


(30) 
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and  £  is  a  4  x  4  matrix  differential  operator.  Equation  (30)  is  cast  as  an 

eigenvalue  equation  by  using  the  functional  form  exp[-ikQSnx],  so 

that  L  e  =  S  e  .  These  e  form  the  set  of  expansion  functions  used  in 
■  ~n  n~n  ~n 

conventional  mode  conversion  theory. 

To  find  an  operator  adjoint  to  a  procedure  discussed  by  Pappert  and 
Smith  (1972)  is  used.  First,  let  the  inner  product  be  defined: 


<a,b> 


b  dz. 


(32) 


where  the  tilde  over  j»  denotes  transpose  and  the  asterisk  denotes  conjuga¬ 
tion.  Next,  the  operator  adjoint  to  is  defined  by  the  relationship 


(33) 


which  is  just  a  prescription  for  integration  by  parts  to  find  the  fora 
of  Jj+  as  well  as  the  necessary  boundary  conditions  on  the  adjoint  functions 
(w)  in  order  that  the  integrated  contributions  vanish.  Pappert  and  Smith 
(1972)  showed  that  the  adjoint  functions  can  be  obtained  from  direct  functions 
by  considering  a  second  waveguide  which  is  related  to  the  original  one  in  a 
simple  way.  In  the  original  guide,  the  z-axis  is  taken  as  vertical  and  propa¬ 
gation  is  assumed  in  the  xz-plane.  The  x-direction-cosine  of  the  geomagnetic 
field  is  l  in  the  original  guide.  The  adjoint  guide  is  obtained  by  replacing 
l  with  -A.  By  examining  symmetry  properties  of  the  elements  of  the 
ionospheric  permittivity  tensor,  Pappert  and  Smith  (1972)  showed  that  the 
adjoint  eigenfunctions  are  given  in  terms  of  eigenfunctions  of  the  adjoint 
waveguide  according  to  the  relationship 
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Further,  from  the  definition  of  adjoint  given  by  equation  (33)  the  usual 
biorthogonality  relations  hold.  Thus 


+ 

<L  W 
■  ~n 


,e  >  -  1*  <w  ,e  > 
~m  n  ~n  "in 

“  <*n*  V 

*  Sm  <W' 


which  shows  both  that  <w  ,e  > 
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0  and  that  X  » 


S*  for  nondegenerate  eigen¬ 


functions. 


MODE  OONVKR8IOH 

Fundamental  to  the  mode  conversion  formulation  employed  here  is  the 
notion  that  any  variation  in  the  VLF  propagation  conditions  along  a  path  can 
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be  modeled  as  a  series  of  discrete  changes  between  homogeneous  sections  along 
the  path.  With  this  in  mind,  it  is  seen  that  the  fundamental  analytical 
requirement  is  the  determi nation  of  mode  conversion  coefficients  due  to  a 
single  abrupt  change  in  waveguide  characteristics. 

Following  classical  procedures,  consider  a  mode  of  unit  amplitude 
incident  on  such  a  discontinuity  in  a  waveguide.  This  mode  is  described  in 
terms  of  its  height-gain  function,  e^fz),  where  the  superscript  I  denotes  the 
incident  region,  the  subscript  j  denotes  the  jth  mode,  the  overarrow  indicates 
a  forward-propagating  mode,  and  the  undertilde  indicates  a  vector  quantity. 
Note  that  the  terms  height-gain  function  and  eigenfunction  are  used  here 
interchangeably  even  though  more  traditional  usage  limits  the  term  height-gain 
to  individual  components  of  the  vector  eigenfunction.  The  discontinuity  in 
the  guide  generates  both  reflected,  or  backward-propagating,  modes  within  the 
incident  region  (I)  and  transmitted,  or  forward-propagating,  modes  within  the 
transmission  region  (II).  Backward-propagating  modes  might  also  exiBt  in 
region  II  due  to  another  discontinuity;  these  will  therefore  be  included.  The 
continuity  requirements  lead  to  the  equation 
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which  is  independent  of  the  position  of  the  discontinuity  between  regions  I 

and  II.  In  equation  (35),  r  .  is  a  coefficient  describing  the  conversion  of 

^  r+I  i 

the  incident  mode  (j)  into  a  backward-traveling  mode  m)  in  region  I  (g^(z)J, 

and  Tni  is  a  coefficient  describing  conversion  of  the  incident  mode  (j)  into 

**  r+Tl  i 

the  forward  traveling  mode  (n)  of  region  II  (z)J.  The  R^,.  are  the  coef¬ 
ficients  appropriate  to  possible  backward-traveling  modes  in  region  II 


Recall  that  a  major  problem  with  implementation  of  mode  conversion 
concepts  arises  from  the  need  to  know  the  future  propagation  history.  In 
terms  of  equation  (35),  this  concept  is  represented  by  the  coefficients  R^. 
By  making  the  assumption  that  appreciable  reflection  does  not  occur,  so  that 
there  can  be  no  backward-traveling  modes,  the  may  be  set  equal 

The  r  .  are  then  also  all  zero,  and  equation  (35)  reduces  to 

*3 


to  zero 


The  assumption  of  negligible  reflection  is  not  unique  to  this  study.  It  is 
commonly  made  and  has  been  demonstrated  to  be  adequate  in  numerous  cases 
(Wait,  1968a,  b;  Galejs,  1971;  Pappert  and  Snyder,  1972}  Smith,  1974).  In 
equation  (36),  there  are  assumed  to  be  N  modes  in  region  II,  so  that  there  are 
N  conversion  coefficients. 

Conventional  mode  conversion  techniques  require  multiplication  of 
equation  (36)  by  the  transpose  conjugate  of  the  adjoint  function  for  the  mth 
forward-traveling  mode  in  region  II,  followed  by  integration  over  all  space. 
Doing  so,  and  recalling  the  biorthogonality  relationships,  we  obtain  the  well- 
known  result 
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more  conveniently  written  as 
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mj  m,m 
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Solving  for  the  conversion  coefficient  yields  the  apparently  simple  result 
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The  sinqilicity  of  equation  (37)  depends  on  knowing  the  height-gain  func¬ 
tions  throughout  all  space.  Such  a  requirement  is  the  main  limitation  to  the 
application  of  mode  conversion  techniques  to  very  long  VLF  paths.  Relaxing 
the  integration  limits  to  some  finite  value  prevents  use  of  the  biorthogonal¬ 
ity  relationships  and  yields,  instead  of  equation  (37),  the  result 
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where  there  are  N  equations  and  N  unknowns,  so  that  the  system  can  be  solved 
for  the  T ny  In  equation  (39), 
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where  the  overarrow  notation  is  dropped  because  only  forward-propagating  modes 
are  considered. 

For  application  to  the  VLF  earth-ionosphere  problem,  the  above  procedure 
may  be  justified  as  follows.  First,  the  various  eigenvalues  and  eigenfunc¬ 
tions  are  determined  by  using  the  coiqplete  geometry  (-•  <  z  <  +“) .  Second, 
the  eigenfunctions  can  be  ejected  to  decay  rapidly  within  the  earth  and 
within  the  ionosphere.  Below  the  ionosphere  but  within  the  waveguide,  the 
vector  eigenfunctions  degenerate  into  the  isotropic  TM  and  TE  components,  with 
a  prescribed  coupling  between  the  two  for  a  given  mode.  By  taking  the  upper 
limit  on  the  integral  to  be  within  the  guide  and  below  the  ionosphere,  it  is 
possible  to  write  the  various  integrals,  1R  in  t  :rms  of  known  functions- 
naroely  the  modififed  Hankel  functions,  as  defined  by  the  Computation  Laboratory 
Staff  at  Cambridge,  Massachusetts  (1945).  The  procedure  used  here  is  actually 
expansion  in  terms  of  "nonorthogonal"  functions  and  is  certainly  better  than 
assuming  orthogonality  and  "throwing  away"  the  contributions  to  the  integrals 
from  the  ionosphere. 

The  jth  vector  eigenfunction  for  either  region  is  defined  to  be  [from 
equation  (31)] 


37 


Because  the  upper  limit  on  the  integrals  in  equation  (40)  is  to  be  taken  below 
the  ionosphere,  equation  (41)  can  be  written  as 
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Because  ^R^(-i)  =  and  jR^(-fc)  =  ^R I  (Budden,  1961b),  equation  (43)  becomes 
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Denoting  this  ratio  as  equation  (42)  becomes 
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The  integrals  in  equation  (44)  are  now  easily  evaluated  by  using  the  con¬ 
ventional  "multiply  and  subtract"  method.  In  this  method,  we  multiply  the 
differential  equation  for  the  height-gain  of  mode  m  by  the  height-gain  of  mode 
n  and  the  differential  equation  for  the  height-gain  of  mode  n  by  the  height- 
gain  of  mode  m  and  subtract.  The  differential  forms  are  obtained  from  equa¬ 
tion  (9a,  b)  and  will  be  presented  shortly  [see  equation  (50)].  Performing  the 
multiplication  and  subtraction  for  hy  components  leads  to 
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The  second  integral  in  equation  (44)  thus  becomes 
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while  the  first  integral  in  equation  (44)  is 
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In  equations  (45)  and  (46),  it  is  assumed  that  C  11  *  C.1  for  all  m  and  j. 

Ul  J 

Evaluation  of  the  second  integral  in  equation  (40)  proceeds  in  the  same 
manner.  The  only  difference  is  the  fact  that  all  field  components  in  this 
case  pertain  to  region  II  of  the  waveguide.  In  this  case  however,  the 
situation  m  *  j  results  in  division  by  zero  in  equations  (45)  and  (46),  and 
special  treatment  is  required.  For  the  case  m  *  j  in  region  II,  the  use  of 
L1 Hospital's  rule  yields 
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(48) 


where  the  superscript  "II”  is  omitted. 

The  contribution  to  the  integral  in  equation  (40)  due  to  the  fields  in 
the  ground  is  obtained  as  follows.  The  fields  vary  as  exp(+ikQTz)  in  the 
ground,  which  gives 
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The  individual  height-gains  are  obtained  in  the  following  way.  Consider 
the  transformed  planar  guide  already  discussed  under  Coordinate  Transforma¬ 
tion.  In  order  to  be  cooqpatible  with  the  assumptions  discussed  under  Mode 
Equation  and  under  Mode  Excitation  Factors  and  Polarization  Coupling,  the 
waveguide  region  is  assumed  to  be  filled  with  the  "equivalent"  permeability  or 
dielectric  medium  as  described  in  under  Coordinate  Transformation.  The  fields 
within  the  “filled*  region  satisfy  equations  (9a,  b)  while  the  fields  within 
the  earth  have  simple  exponential  forms. 

Since  all  altitudes  of  interest  are  small  compared  to  the  radius  of  the 

earth  (a),  the  exponential  term  in  equation  (9)  is  expanded  to  first  order  in 

~i.k  Sx 

z/a.  A  variation  in  the  x-direction  of  the  form  e  °  is  assumed.  Note 
that  the  variable  S  will  be  identified  later  as  sin(6),  where  the  angle  6  is 
referenced  to  the  vertical  rather  than  to  the  horizontal  as  in  previous 
sections.  With  the  above  assumptions,  equation  (9b)  for  TM  waves  becomes 


Hy  +  (C2  +  az)  Hy  =  0, 


(50) 


where  C2  *  1-s2,  a  -  2/a,  and  the  double  prime  denotes  second  differentiation 
with  respect  to  z. 

With  the  change  of  variable 
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the  equation  becomes 


H"  +  tH  *  0,  (52) 
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where  now  the  double  prime  denotes  differentiation  with  respect  to  t. 
Equation  (52)  is  Stokes'  equation.  Thus,  is  linearly  expressible  in  terms 
of  solutions  to  Stokes'  equation.  In  this  study,  we  shall  use  the  modified 
Hankel  functions  of  order  one-third  ( h 1 , h2 ) ,  as  defined  by  the  Computation 
Laboratory  Staff  at  Cambridge,  Massachusetts  ( 1945) .  Using  a  similar  pro¬ 
cedure  for  TE  polarization,  we  shall  arrive  at  similar  results  except  that  the 

transformed  planar  guide  will  contain  a  medium  with  refractive  index  given  by 
ez/a. 

The  waves  in  the  two  regions  are  given  by  the  following  relationships: 


Region  I: 
z  >  d 


Hy  =  a^hj(t)  +  a^h^tt) 

*«-  l[£fr3 


Region  II: 
z  >  0 


H 


y 


ik  Tz 
o 
e 


E 

X 


T 


ik  tz 
o 

e  , 


42 


whore  hj  and  are  modified  Hankel  functions  of  order  one-third,  N  is  the 
refractive  index  of  the  ground,  and  T  -  -  S  J  '  ,  taken  with  negative 

imaginary  part.  Continuity  requirements  between  the  regions  leads  to  the 
results 
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Following  a  similar  procedure,  results  for  the  TE  polarization  are  given: 
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Normalizing  the  Hy  component  to  1  at  z  =  0  and  using  equation  (28)  for  the 
polarization  coupling,  we  arrive  at  the  final  forms  for  the  height  gains: 
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Once  the  conversion  coefficients  are  obtained,  it  is  necessary  to  sum 
over  all  the  incident  modes.  This  requires  evaluation  of  the  amplitudes  of 
the  incident  modes.  This  can  be  done  by  starting  at  the  transmitter,  where 
the  modal  excitation  factors  give  the  mode  amplitudes,  and  by  assuming  that 
the  regions  between  the  discontinuities  are  homogeneous. 

Notations lly  this  can  be  written  as 
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where  the  equation 
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represents  the  angilitudes  of  the  incident  nodes  (a*]  at  the  end  of  the  region 
I  and  the  equation 
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represents  the  amplitudes  of  the  modes  at  the  start  of  region  I.  The 
factor  g1  is  a  diagonal  matrix  with  the  nth  diagonal  element  given  as  follows: 
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Note  that  Pn  represents  the  propagation  factor  appropriate  to  a  homogeneous 
section  of  length  t  Ap ) 1 ,  assumed  to  be  the  length  of  region  1.  Note  that  the 
amplitudes  of  the  modes  in  region  II  just  after  conversion  can  be  conveniently 
written  as 
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(58) 


where  J  is  a  rectangular  matrix  of  the  conversion  coefficients  from  region  I 
into  region  II . 

A  mode  conversion  formulation  has  been  presented  by  Pappert  and  Shockey 
(1974),  but  important  differences  exist  between  their  treatment  and  that 
engtloyed  here.  Differences  occur  both  in  the  numerical  procedures  enployed 
and  in  the  definition  of  the  conversion  coefficients.  In  the  present  formula¬ 
tion  a  conversion  coefficient  is  computed  at  each  slab  interface  independently 
of  any  previous  conditions.  In  the  Pappert  formulation,  the  conversion  coef¬ 
ficients  are  defined  and  computed  in  such  a  way  that  the  effect  of  the  propa- 
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gat ion  factor  [equation  (57))  is  included.  Both  definitions  of  conversion 
coefficients  are  reasonable,  although  their  maaerical  values  can  be  quite 
different.  For  example,  first  consider  an  isolated  terminator  specified  as  an 
abrupt  transition.  The  two  formulations  will  give  the  same  results  for  the 
conversion  coefficients.  Next  suppose  that  the  terminator  has  a  substantial 
extent.  The  conversion  coefficients  then  differ  for  each  slab  interface  after 
the  first  because  the  Pappert  formulation  includes  the  accumulated  effects  of 
propagation  through  the  transition  region,  whereas  the  current  formulation 
does  not.  In  any  case,  the  resulting  mode  sums  computed  by  means  of  either 
formulation  will  be  sensibly  the  same,  with  differences  primarily  due  to 
confute r  numerics. 

The  Pappert  definition  of  conversion  coefficients  is  useful  for  studying 
a  horizontal  inhomogeneity  that  exists  over  a  small  distance,  such  as  the 
dawn-dusk  terminator.  A  similar  conversion  coefficient  can  be  derived  from 
the  mode  constants  and  the  mode  conversion  coefficients  presented  in  this 
study.  What  is  required  is  the  complex  amplitude  of  a  mode  at  the  end  of  such 
a  confined  inhomogeneity  resulting  from  a  single  mode  of  unit  amplitude  at  the 
start  of  the  inhomogeneity. 

Denote  by  region  1  the  region  prior  to  the  start  of  the  inhomogeneity. 
The  complex  amplitudes  of  the  modes  in  the  first  region  of  inhomogeneity 
(called  region  2)  due  to  a  mode  of  unit  amplitude  and  number  k  in  region  1  are 

2  i 

given  by  the  conversion  coefficients  T  '  [see  equation  (44)] .  The  complex 

m,  k  2 

amplitudes  of  the  modes  at  the  end  of  region  2  are  A  ,  where 

m,  k 


2 

m,k 


2 

and  PB  is  the  propagation  factor  for  mode  m  in  region  2  [see  equation  (57)]. 

The  complex  amplitude  of  mode  n  at  the  start  of  the  next  region  (region  3)  due 

3  2 

to  a  mode  of  unit  amplitude  in  region  2  is  T  '  ,  so  that  the  actual  complex 

n,B  3  2 

amplitude  of  the  mode  n  in  region  3  is  given  by  the  product  of  T  '  end 

n,m 

.  Summation  over  all  the  modee  of  region  2  gives  the  total  complex 
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amplitude  of  mode  n.  Thus 
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(59) 


This  result  can  be  generalized  to  any  inhomogeneity  along  the  path.  We  obtain 
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The  formulations  presented  in  this  section  have  been  implemented  in  a 
computer  program  written  in  FORTRAN.  The  development  of  the  numerical  cal¬ 
culations  required  for  the  reflection  coefficients  and  waveguide  modes  has 
occurred  over  many  years.  These  efforts  have  been  supported  by  the  Defense 
Nuclear  Agency,  NAVEL  EX  PME  117,  the  Defense  Communications  Agency,  and  the  OS 
Coast  Guard. 
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4.  name 


FAST MC  is  a  computer  program  which  employs  the  approximate  mode  conver¬ 
sion  model  described  in  section  3.  The  program  is  written  in  Univac's  version 
of  ASCII  FORTRAN.  It  accepts  card  image  input  and  generates  plots.  Addition¬ 
al  outputs  available  are  printouts  of  the  amplitude  and  phase  and  binary  out¬ 
put  of  the  complex  field.  The  source  is  an  arbitrarily  oriented  dipole  at  any 

altitude  within  the  waveguide.  The  vertical  (e_)  and  horizontal  broadside 

z 

(ey)  fields  at  any  altitude  can  be  obtained. 

Input  to  the  program  is  controlled  by  a  set  of  control  cards  which  define 
the  type  of  data  being  read  and  allow  standard  defaults  to  apply.  These  con¬ 
trol  cards  are  described  below.  All  eighty  columns  of  these  cards  are  read 
and  printed.  However,  only  the  first  four  columns  are  examined  by  the  pro¬ 
gram.  Thus,  the  input  card  and  the  printout  can  contain  additional  comment. 
For  example,  "NAMELIST  DATA  FOR  ELEVATED  TRANSMITTER"  can  be  used  in  place  of 
the  minimum  requirement  of  "NAME".  The  control  cards  are  described  below,  and 
a  sample  set  is  shown  in  figure  1. 

NAME  Signals  that  NAMELIST  data  follow.  In  this  program  the 

NAMELIST  name  is  DATUM.  The  NAMELIST  variables  are 
described  below. 

I COMP  =  Index  of  received  electric  field  conponent 

=  1  gives  the  vertical  field  (Z)  -  DEFAULT 
-  2  gives  the  horizontal  transverse  field  (Y) 

TALT  *  Transmitter  altitude  in  Ion  -  DEFAULT  *  0.0 

RALT  -  Receiver  altitude  in  km  -  DEFAULT  *  0.0 

INCL  »  Inclination  of  transmitting  antenna  with  respect  to  the 

vertical  in  degrees  -  DEFAULT  =*  0.0. 

THETA  «  Orientation  of  transmitting  antenna's  projection  in  the 

horizontal  plane  with  respect  to  the  direction  of  propa¬ 
gation  in  degrees  -  DEFAULT  m  0.0. 
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NAME 

40ATUM 

NPRINT*3. 

DMAX*  3 . ,  SIZEX1  *6.  , 

SEND 

DATA 

FASTMC  TEST,  0-0 

R  0.0  f  31.1350  A  208.334  C  4.3:18  M  3-600D-05  S  4.640D  00  E  81.0 

1  89.78599  -6 . 2589  22-8.0 1 75551 90-06-3 . 14370439D-06-1 . 75676358D-1 2-6 . 2379691  ID-1  3 

2  89.78599  -6 . 258922-3. 6C885085D-09-1 . 8062 H I 8D-09  9. 48673970D-01 -1 . 861 1 4244D-0 1 

1  88 . 984  t  3  -1.682941  4 . 552990030-03- 1 . 98787729D-02-2 . 7 1 479704D- 1 2  2 .6551 1 284D-1 2 

2  88.98413  -1.682941-2.440345930-07-1.551706140-07  9.594997620-01-1.699363550-01 

1  83.31784  -0.608891  1 . 25067533D-03-2 . 3366 1 3450-02-3 . 8 1 5659600- 1 2  2.185022450-12 

2  83.31784  -0 . 608891 -2 . 60604588D-07- 1 . 909885790-07  9. 7 1 503 1 32D-01 - 1 . 463 1 55000-0 1 

1  78.98910  -0.665431  9 . 73001 6340-04-?. 02983 1 050-02-1 . 55572876D-1 1 -3 . 427403330-1 2 

2  78.98910  -0.665431-3.2205333‘-'0-07-4.715668900-07  9.901953500-01-1.042478400-01 

1  74.97312  -0.764111  1 . 234B7260D-04- 1 . 967820050-02-1 . 086145450-1 1 -2 .846834350-1 1 

2  74.97312  -0 . 7641 1 1 - 1 . 1 9 1 89454D-07-7 . 552847030-07  1.010702380  00-3.736264340-02 

R  1.000  F  31.1350  A  208.234  C  4.308  M  5.600D-05  S  1.0000-05  E  5.0 

1  89.68035  -5.566901-1 .24182979D-05  1.278074990-04  7.365554290-09  2.507921680-08 

2  89.68035  -5.566901  4 . 1 42 1 029 1 D-09- 1 . 856645390-06  9. 5073475 1 D-01-1 . 821 507030-0 1 

1  86.05721  -1.194791  4 . 4863599' 0-04  1.682497450-03  7.055075460-08  2.618027010-07 

2  86.05721  -1 . 1947U1-6.50431216D-06-2. 1472t 1960-05  9.616847730-01-1.601697970-01 

1  80.93273  -1.269851  2.843831290-03  5 . 42606425D-03  2.449055440-07  5.717887800-07 

2  80.93273  -1 . 269851 -2 . 27 1 62 1 20D-05-5 . 700980450-05  9. 7464841 30-01 -1 . 23476674D-01 

1  76.61480  -1.636121  1.082251110-02  1.109913520-02  1.935805520-08  8.073643290-07 

2  76.61480  -1 .636121-4.015206300-05-1 . 05085694D-04  9 . 87201 1 64D-01 -6 . 430890030-02 

1  72.32836  -2.025471  3.268147330-02  6 . 0003 1 0250-03-5 . 4601 565 1 D-07  7 . 30226214D-07 

2  72.32836  -2.025471-6.720166230-05-1.616247640-04  9 .87 1 700B3D-01  2 .67031 249D-02 

R  1.500  F  31.1350  A  208.234  C  4.398  M  5-600D-05  S  4.6400  00  E  81.0 

1  89.78599  -6 . 258922-8.0 1 7555190-06-3 . 1 43704390-06-1 . 756763580-12-6 . 2379691 1 D-1 3 

2  89.78599  -6.258922-3.668850850-09-1 .80621 1 18D-09  9 . 48673970D-01-1 . 861 1 42440-01 

1  88.98413  -1.682941  4 . 55299003D-03- 1 . 9B787729D-02-2 . 7 1 4797040-1 2  2.655112840-12 

2  88.98413  -1.682941-2.440345930-07-1.551706140-07  9 . 59499762D-01 -1 . 69936355D-01 

1  83.31784  -0.608891  1.250675330-03-2.336613450-02-3.815659600-12  2.185022450-12 

2  83.31784  -0.608891-2.606045880-07-1.909885790-07  9.715031320-01-1.463155000-01 

1  78.98910  -0.665431  9-730016340-04-2.029831050-02-1.555728760-11-3.427403330-12 

2  78.98910  -0.665431-3.220533350-07-4.715668900-07  9.90195350D-01-1 . 04247840D-0 1 

1  74.97312  -0.764111  1 . 234872600-04-1 . 967820050-02-1 . 086145450-1 1-2 .846834350-1 1 

2  74.97312  -0.764111-1.191894540-07-7.552847030-07  1.010702380  00-3.736264340-02 

R  40. 

NAME 

SDATUM  NPRINT-1,  TALT-6.,  RAlT-6.,  INCl-90.,  THETA-90.,  4EN0 
START 
NAME 

40ATUM  INCL«0. ,  THETA-O.,  NPPLOT-1,  NPRINT-O,  1C0MP-2,  .END 
START 
NAME 

40ATUM  INCL-90.,  THETA-90.,  «END 
START 


Figure  1.  Sample  card  deck  for  FASTMC. 
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TOPHT 

NPRINT 


NAP LOT 

AMP MIN 

AMPMAX 

AMPINC 

NPPLOT 

PHSMIN 

PHSMAX 

PHSINC 


*  Height  of  the  bottom  of  the  ionosphere  in  Ian  [integration 
limit  h  in  equations  (49)  -  (53))  -  DEFAULT  *  90.0 

“  Flag  to  control  printout  of  mode  conversion  and  mode  sum 
results 

*  0  gives  no  printout 

**  1  gives  printout  of  mode  parameters  and  calculated  fields 
vs.  distance  -  DEFAULT 

*  2  also  prints  the  input  mode  constants  cards  and  the  mode 
conversion  coefficients 

*  3  also  prints  the  integrals  used  in  the  calculation  of 
the  mode  conversion  coefficients. 

=  Flag  to  control  amplitude  plots 

31  0  gives  no  plot  of  aag>litude 

=  1  gives  amplitude  plots  -  DEFAULT 

*  Minimum  amplitude  to  be  plotted  in  dB  above  one 
microvolt/meter  -  DEFAULT  =  10.0. 

«  Maximum  anqalitude  to  be  plotted  -  DEFAULT  =  70.0 

=  Tic  mark  interval  (for  plot  border)  -  DEFAULT  =  10.0 

-  Flag  to  control  plots  of  phase  relative  to  the  speed  of 
light 

*  0  gives  no  plot  -  DEFAULT 

»  1  gives  phase  plots 

=  Minimum  phase  excursion  in  degrees  -  DEFAULT  *  360.0 

x  Maximum  phase  excursion  -  DEFAULT  -  360.0 

-  Tic  mark  interval  (for  plot  border)  -  DEFAULT  ■  90.0 


50 


DMIN 


=  Minimum  distance  from  the  transmitter  in  Mm  -  DEFAULT 


0.0 

DMAX  =  Maximum  distance  -  DEFAULT  =*  10.0 

XINC  =  Horizontal  tic  mark  intervals  (for  plot  border)  -  DEFAULT 

-  1.0 

Length  of  horizontal  axis  in  inches  -  DEFAULT  =  10.0 

Length  of  vertical  axis.  Maximum  value  is  9.0  -  DEFAULT 

=  6.0 

Number  of  sets  of  data  to  be  plotted  in  one  graph. 
Maximum  value  is  4  -  DEFAULT  =  1 

Number  of  data  points  to  be  calculated.  Maximum  value  is 
501  -  DEFAULT  *  501 

TOT APE  »  Binary  output  flag 

«*  0  gives  no  output  -  DEFAULT 

=  1  writes  the  complex  mode  sum  and  path  parameters  to 
logical  unit  2 

TLONG,TLAT  =  Transmitter  coordinates  in  degrees  west  and  north  - 
DEFAULT  *  0.0,  0.0 

RBEAR  *  Path  bearing  in  degrees  clockwise  from  north  -  DEFAULT  = 

0.0 

After  reading  the  NAMELIST  the  program  reads  the  next  control  card. 

DATA  Signals  that  the  propagation  path  data  follow.  The  format 

of  this  data  is  that  which  is  produced  by  MDDESRCH  (Morfitt 
and  She liman,  1976),  WAVEGUID  (Pappert  et  al,  1970)  and  the 
INTEGRATED  PREDICTION  PROGRAM  (Ferguson,  1972).  The  first 


SIZEX 

SIZEY 

NRCUKV 

NRPTS 
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card  contains  the  data  set  identification.  This  card  is 
used  to  supply  a  literal  path  description  such  as  "HAWAII 
TO  WAKE".  This  is  followed  by  sets  of  mode  constants,  one 
for  each  path  segment  or  slab.  The  first  card  for  each 
segment  contains  the  starting  distance  (P),  frequency  (f), 
magnetic  azimuth,  codip  and  intensity,  ground  conductivity 
and  permittivity  (0,  £),  and  the  nominal  height  of  the  free 
space  portion  of  the  waveguide  (h).  Of  these  parameters, 
P,  f,  o,  and  e  must  appear.  The  magnetic  parameters  are 
included  for  reference.  The  value  of  h  may  be  included  if 
it  is  to  vary  among  the  segments.  Otherwise  it  is  a  con¬ 
stant  value  input  via  the  NAMELIST.  This  card  is  followed 
by  the  mode  constants  cards,  one  pair  for  each  mode,  con¬ 
taining  the  following  information: 


1  6  I  T,  t2 

2  6  i  T3  T4 

where  1  and  2  are  sequencing  indices,  6  is  the  complex 
eigenangle  at  the  ground,  I  is  an  index  indicating  the 
functional  form  of  the  ratio  ey/hy  [equation  (28)]  and  the 
Ts  are  complex  constants  described  below. 

The  list  of  modes  for  each  segment  is  terminated  by  a  blank 
card.  A  maximum  of  20  modes  will  be  used  although  there  is 
no  maximum  number  which  the  data  deck  nay  contain.  The 
list  of  segments  is  terminated  by  a  card  with  P  =  40. 
Parameters  for  each  segment  are  stored  sequentially  on 
logical  file  4.  The  number  of  path  segments  is  limited  by 
the  space  allocated  to  this  file.  Each  segment  requires 
1848  words  of  storage.  After  reading  these  data  the  program 
generates  the  required  calculations  (as  specified  by  the 
NAMELIST  input  or  defaults)  then  returns  for  another  con¬ 
trol  card. 
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BUMP 


Signals  that  this  card  contains  an  overall  plot  iden¬ 
tification.  It  is  used  for  additional  information  not 
contained  in  the  data  set  identification  cards.  When 
several  sets  of  data  are  plotted  on  one  graph,  it  can  be 
used  to  provide  a  connection  between  them. 

START  Signals  calculations  using  new  NAMELIST  parameter  values 
for  the  current  propagation  path  data. 

A  summary  of  how  the  control  cards  are  used  in  FASTMC  is  shown  in  figure 
2.  A  brief  description  of  the  subroutine  functions  is  given  in  table  1,  and 
complete  listings  can  be  found  in  appendix  A.  Sample  plotted  outputs  are 
shown  in  figures  3  and  4.  The  Ts  on  the  data  cards  are  given  by  the  following 
relationships: 
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MCFIDS:  Calculate  and  store  mode  sum  vs  distance 


MCPLTS:  Print/plot  mode  sum 


i) 

Figure  2.  Illustration  of  FASTMC  program  execution  sequence- 
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MCSTEP  Calculates  the  integrals  (50)  and  solves  for  the  conversion 

coefficients  (44)  at  each  segment  boundary.  Calls  MDHNKL 
and  CLINEQ. 

MCFLDS  Calculates  the  field  strength  from  DMIN  to  DMAX.  Calls 

HTGAIN. 

MCPLTS  Prints  and  plots  amplitude  and  phase.  Calls  BORDER,  C3JRVE, 

and  CALCOMP  routines. 

MDHNKL  Calculates  modified  Hankel  functions  of  order  one-third  and 

their  derivatives. 

HTGAIN  Calculates  height-gain  functions  (32,  33).  Calls  M)HNKL. 

CLINEQ  Solves  system  of  equations  for  complex  coefficients. 

BORDER  Draws  border  for  plots.  Labels  minimum  and  maximum  values 

of  axes. 

CURVE  Draws  a  curve  of  x  vs  y.  Seven  different  line  types  are 

available.  FASTMC  uses  only  the  first  four  types:  solid, 
long,  medium,  and  short  dashes. 


Table  1.  FASTMC  subroutines. 


Z  COMPONENT  AMPLITUDE  POWER  =  1 .  FREQ  =  31.135 
I NCL  =  0.  THETA  =  0.  TALT  =  0-0  RALT  =  0-0 

FASTMC  TEST.  0=0 


Figure  3.  Sample  plot  of  amplitude  output  from  FASTMC.  This  is 
the  first  plot  generated  by  the  card  deck  of  figure  1 


Prom  these  relationships,  equation  (26)  can  be  written  as  follows: 
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and  equation  (28)  can  be  evaluated  at  the  ground,  giving  the  relationship 
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where  f(,  f ,  are  from  equations  (53)  and  (54). 

The  relationship  in  equation  (62)  is  used  if  the  value  of  I  on  the  mode 
cards  is  1 .  If  the  value  of  I  is  2,  the  following  relationship  is  used: 


Y  = 


(i  +  jijKi  -  |R,  ,R,) 
(i  +  ,R|)  L*L  j.R| 


f,(d) 

f][TdT 


T  T 
3  4 


(63) 


Equations  (62)  and  (63)  are  identical  if  the  value  of  the  mode  equation  (12) 
is  exactly  zero.  The  value  of  I  is  determined  from  the  mode  polarization  as 
calculated  in  the  WAVBGUID,  (CDKSRCH ,  or  I NT EG RATH)  PREDICTION  PROGRAM.  A 
mode  which  is  principally  TM  will  have  1  -  ^R^  ^R^  very  small  and  the  ratio  Y 
will  be  very  large,  resulting  in  the  use  of  equation  (62).  Conversely,  a  mode 
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which  is  principally  TE  will  have  1  -  (R(  (R f  very  small,  and  the  ratio  Y  will 
be  very  small,  resulting  in  the  use  of  equation  (63). 

The  excitation  factor  for  due  to  a  vertical  radiator  can  be  defined 
from  equation  (61)  as  -ST^.  Excitation  factors  due  to  horisontal  end-fire  and 
horizontal  broadside  radiators  can  be  obtained  from  Pappert  et  al,  1970.  In 
terms  of  the  Ts ,  these  are  T1  and  T3T4  for  end-fire  and  broadside,  respective¬ 
ly.  The  corresponding  electric  field  conponents  are  obtained  from  equation 
(55). 
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BUM>  is  a  computer  program  which  has  been  written  to  allow  for  a  system¬ 
atic  variation  of  a  path  perturbation,  with  full  allowance  for  variation  of 
the  mode  parameters  with  conductivity  and  magnetic  field.  The  most  common 
perturbation  is  the  day-night  terminator.  Previous  studies  of  the  effect  of 
the  terminator  on  propagation  were  limited  to  waveguides  for  which  the  mag¬ 
netic  field  and  ground  conductivity  were  constant.  BUMP  allows  for  the  study 
of  the  effect  of  the  terminator  along  a  path  with  mixed  ground  conductivities 
and/or  significant  magnetic  field  variations  along  the  propagation  path.  The 
formulation  of  the  program  allows  for  other  perturbations  to  be  studied. 

The  design  of  the  program  was  aimed  at  requiring  the  simplest  possible 
data  deck  preparation.  This  was  done  in  order  to  make  it  easy  to  modify 
propagation  path  data  and  to  allow  flexibility  in  the  ionospheric  models.  The 
resulting  input  requires  propagation  path  data  for  a  series  of  horizontally 
homogeneous  ionospheres.  The  data  for  each  ionosphere  are  merged  by  the 
program  to  produce  the  required  ionospheric  perturbation.  This  program  has 
been  used  to  study  models  of  the  ionospheric  variation  from  middle  to  high 
latitudes  and  to  study  trans-terminator  propagation  across  the  magnetic 
equator.  The  program  to  be  described  below  is  set  up  for  a  minimum  path 
segment  length  of  5  km  and  a  maximum  path  length  of  40  Mm.  This  requires  801 
random  access  records  per  ionosphere.  Since  a  substantial  portion  of  the  cost 
of  running  the  program  is  incurred  in  setting  up  the  random  access  files 
required,  users  can  achieve  some  cost  savings  by  changing  the  5  km  and  40  Mn 
figures  and  modifying  the  random  access  files  accordingly.  A  maximum  of  30 
modes  can  be  stored  for  each  path  segment,  requiring  1229  words  of  storage  in 
each  random  access  record. 

Preparation  of  the  input  data  for  this  program  is  illustrated  for  a  bump 
in  the  ionosphere.  The  lowest  height  of  ionosphere  will  be  characterized  by 
h*  ■  hQ  and  the  highest  height  by  h'  »  hg.  The  variation  of  the  ionosphere 
from  lowest  to  highest  is  modeled  by  a  linear  variation  of  h'  over  six  steps 
denoted  by  h^,  h2»  h3,  h4 ,  hg.  The  length  of  each  step  is  Ap^  and  the  length 
of  the  step  at  h*  is  Ap2*  The  bu^»  in  the  propagation  path  is  defined  in 
figure  5. 


60 


Illustration 


By  making  Ap.  equal  to  20  Mm,  the  bump  becomes  a  ramp  for  modeling  the 
day-night  transition.  If  the  order  of  input  of  the  data  sets  corresponding  to 
each  h'  were  reversed,  the  ramp  would  model  the  night-day  transition.  Movement 
of  the  bump  is  accomplished  by  varying  p  from  pQ  to  PL  in  increments  of  Ap. 
If  Pq  is  greater  in  value  than  PL,  the  sign  of  Ap  is  adjusted  automatically. 
The  total  path  length  is  P^j,* 

If  the  data  set  for  hR  does  not  have  data  corresponding  to  P  +  (n  -  1) 
Ap.j,  the  data  from  the  next  smaller  distance  is  used.  For  example,  suppose  hR 
has  data  at  1.0  Mn  and  1.2  Ita  and  p  ♦  (n  -  1)Ap^  i«  greater  than  1.0  Mn  and 
less  than  1.2  Mm.  The  program  will  select  the  data  for  1.0  Mm  to  be  used  at  P 
+  (n  -  1)Ap  .  if  the  value  of  p  becomes  less  than  zero,  the  necessary  steps 
in  hn  are  removed.  For  example,  if  p  =  -1  •5Ap1 ,  the  first  step  on  the  propa¬ 
gation  path  will  be  h^.  A  similar  procedure  is  followed  if  p  +  (n  -  1)Ap^ 
becomes  greater  than  pMX. 

The  input  to  BUMP  is  a  series  of  FASTMC  compatible  data  sets,  one  for 
each  ionospheric  profile.  Each  data  set  must  contain  data  for  the  required 
propagation  path.  The  output  from  BUMP  is  a  series  of  FASTMC  data  sets,  each 
corresponding  to  one  position  of  the  bump  or  ramp  on  the  propagation  path. 
The  data  set  identification  card  which  follows  each  "DATA"  card  contains  the 
values  of  p,  Ap1 ,  Ap2,  and  the  hs  heights.  An  additional  identification  card 
corresponding  to  the  FASTMC  input  control  card  "BUMP*  may  be  output. 

The  order  of  execution  of  the  program  is  specified  by  a  set  of  control 
cards.  Each  card  is  read  and  printed  in  its  entirety,  but  only  the  first  four 
columns  are  examined  by  the  program.  These  control  cards  and  their  functions 
are  described  below. 

NAME  signals  that  NAMELIST  data  follows.  The  NAMELIST  is  DATUM.  The 
input  variables  are  defined  as  follows: 


RHOO 


Pp  in  (ta  -  DEFAULT  =  1.0 


RHOL  =  PL  in  Mm  -  DEFAULT  =2.0 
DRHO  =  Ap  in  fta  -  DEFAULT  =  0.1 
DRHOI  =  Apt  in  m  -  DEFAULT  =  0.1 
DRH02  =  Ap2  in  fti  -  DEFAULT  =  20.0 
RHOMAX  =  pffiax  in  Mn  -  DEFAULT  =  10.0 

HPRIME  =  h'  in  km.  Maximum  number  of  values  is  7  -  DEFAULT  = 
87.0,  85.0,  83.0,  81.0,  79.0,  77.0,  75.0 

ID  signals  that  the  next  card  is  to  be  read  for  identification  of 
the  output.  The  first  76  columns  of  this  card  are  stored  on  a 
card  with  BUMP  in  its  first  4  columns.  This  can  be  processed 
by  FASTMC . 

DATA  signals  that  data  corresponding  to  one  of  the  values  of  h' 
follow.  These  data  are  in  the  same  format  as  required  by 

FASTMC.  They  oust  be  input  in  the  order  indicated  by  the 

HPRIMEs  in  the  NAMELIST.  The  data  for  HPRIME (n)  are  stored  in 
logical  file  10  +  u.  The  individual  segments  are  stored  in 
random  access  files,  with  each  record  representing  a  multiple 
of  5  km  with  a  maximum  distance  of  40  Mm.  Thus,  data  at  1.234 
Mm  will  be  stored  as  1.235  Mm  in  record  number  248. 

START  signals  that  all  necessary  inputs  are  coaqplete  and  the  genera¬ 
tion  of  the  output  data  sets  is  to  begin.  The  output  is 
written  to  external  file  "OUTPUT"  (logical  file  2)  in  FASTMC 
card  image  format. 

The  resulting  data  are  used  in  FASTMC  to  obtain  the  variation  of  field 
strength  parametric  in  the  position  of  the  perturbation.  Various  programs 
exist  or  could  be  written  to  display  this  output  in  other  forms.  Appendix  B  is 
a  listing  of  the  program. 
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To  illustrate  the  execution  of  BUMP  with  FASTMC,  let  us  assume  that  a 
file  named  NLK  exists*  This  file  contains  data  elements  for  a  variety  of 
ionospheric  profiles  on  the  propagation  path  NLK  to  Brisbane,  Australia*  The 
electron  density  profile  is  exponential  in  height  with  a  slope  B  and  a  refer¬ 
ence  height  h'  denoted  by  8/h*.  The  data  elements  in  NLK  are  renamed  BBHH, 
where  BB  and  HH  denote  8  and  h' ,  respectively*  For  exaiqple,  the  daytime 
ionosphere  is  taken  to  be  0.3/74  and  the  corresponding  data  element  is  named 
3074*  The  nighttime  ionosphere  is  taken  to  be  0.5/86.  The  variation  in 
ionosphere  from  day  to  nicflit  is  made  in  five  equally  spaced  steps  in  B  and 
h' .  The  output  from  BUMP  is  to  go  into  a  temporary  file  named  OUTPUT,  which 
is  to  be  used  in  FASTMC  to  generate  the  resulting  mode  sums*  The  data  deck  to 
carry  out  these  operations  is  shown  in  figure  6.  The  resulting  plot  is  shown 
in  figure  7.  Note  that  the  use  of  the  ”10"  control  card  in  the  input  to  the 
BUMP  program  results  in  additional  identification  of  the  curves  via  the  "BUMP” 
control  card  in  FASTMC.  Also  note  that  although  B  as  well  as  h*  varies,  there 
is  provision  for  identifying  only  the  latter;  the  B  variation  must  be  indica¬ 
ted  in  the  "ID"  card. 
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SASG.T  OUTPUT ., F/// 1000 

exgr  p  ,v.  bump 

id 

NLK  TO  BRISBANE,  NIGHT  ( 0 . 5/86  >  T0  DAY  (0.3/74) 

NA'.'E 

SDATUU  RHOO=2.,  RH0L*5.,  DRH0  = ' • •  DRH01=.150, 

R  HO-.' A'.  =  12.  , 

HPR I ME = 86 . ,  84 . ,  82.,  SO..  7B ■ .  76.,  74., 

SEND 

PADD  NLK.508G 
3»ADD  NLK.  4784 
SADD  NLK. 4382 
S»ADD  NLK. 4080 
PADD  NLK. 3778 
SADD  NLK. 3376 
@A00  NLK. 3074 
START 

0XQT  PGM . FASTMC 
NAME 

SDATUM  NPRINTaO,  NRCURV.4,  D*|AX*12.,  SIZEX=6.,  AMPMIN  =  -10.,  SIZEY=4.,  SEND 
•ADO  OUTPUT. 


Figure  6,  Sample  card  deck  for  running  BUMP  and  FASTMC. 
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subroutine  TO  CALCULATE  MODE  conversion  coefficients 
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